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I.  Statement  of  Work 


Many  United  States  Air  Force  routing  and  scheduling  problems  can  be  analyzed  and 
modelled  mathematically  as  some  type  of  network  flow  problem.  As  budget  pressures 
become  more  intense,  there  will  be  even  more  emphasis  on  optimizing  the  use  of  the 
assets  available  to  the  Air  Force.  Optimization  algorithms  and  state-of-the-art  software 
will  be  important  tools  that  can  be  used  by  Air  Force  planners  and  decision-makers 
during  the  foreseeable  future.  Models  and  techniques  which  exploit  the  underlying 
network  structure  of  these  problems  will  also  gain  added  importance.  The  research 
reported  in  this  document  is  directed  toward  the  discovery  of  new  and  improved 
algorithms  and  the  development  of  computational  software  to  assist  in  the  area  of  routing 
and  scheduling  problems. 


II.  Technical  Reports 

Title 

Network  Flows 
(September  1992) 

Authors 

Richard  V.  Helgason 
Jeffery  L.  Kennington 

Executive  Sumary 

The  objective  of  this  report  is  to  clearly  present  the  techniques  used  in  a 
computationally  efficient  implementation  of  the  primal  simplex  algorithm  specialized 
for  various  network  models.  Numerous  examples  are  included  to  help  illustrate  the 
ideas  and  techniques.  It  is  shown  mathematically  and  illustrated  by  example  that 
every  basis  for  a  linear  network  problem  is  triangularizable  and  corresponds  to  a 
spanning  tree.  This  result  is  used  in  the  development  of  a  specialization  of  the  primal 
simplex  algorithm  for  the  linear  network  problem  which  is  at  least  two  orders  of 
magnitude  faster  then  general  purpose  software.  It  is  shown  mathematically  and 
illustrated  by  examples  that  every  basis  for  a  generalized  network  problem  has  a  block 
diagonal  structure  with  each  block  corresponding  to  either  a  rooted  tree  or  a  one-tree. 
While  this  structure  is  not  as  beneficial  as  that  of  the  linear  network  problem,  it  can 
still  be  exploited  in  performing  the  simplex  operations.  Other  similar  results  for  the 
multicommodity  network  flow  problem  are  also  presented. 

Publication  Status 

This  manuscript  has  been  submitted  for  publication  as  a  chapter  in  a  book  and  is 
currently  under  review. 
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Title 

The  One-to-One  Shortest-Path  Problem:  An  Empirical  Analysis  with  the  Two-Tree 
Dijkstra  Algorithm 
(Revised  December  1992) 

Authors 

Richard  V.  Helgason 
Jeffery  L.  Kennington 
B.  Douglas  Stewart 

Executive  Summary 

The  problem  of  finding  the  shortest  path  between  a  designated  pair  of  nodes  in  a 
graph  is  a  fundamental  problem  in  operations  research,  computer  science,  and 
scheduling  theory.  Good  algorithms  for  this  problem  can  be  used  as  a  building  block 
for  other  algorithms  to  solve  the  assignment  problem,  the  semi-assignment  problem, 
the  multicommodity  network  flow  problem,  and  integer  networks  with  side  constraints. 
The  classical  Dijkstra  algorithm  begins  at  one  of  the  designated  nodes  and  fans  out 
from  this  node  until  the  other  designated  node  becomes  a  member  of  the  labeled  set. 
This  is  easily  accomplished  by  a  computer  implementation  that  builds  a  tree  rooted  at 
one  of  the  designated  nodes.  In  our  investigation  we  empirically  demonstrate  that  a 
better  algorithm  (in  terms  of  computational  time)  is  obtained  by  a  procedure  that 
begins  at  both  designated  nodes  and  fans  out  in  both  directions,  either  simultaneously 
as  in  our  parallel  implementation  or  alternately  as  in  our  sequential  implementation. 
This  new  algorithm  terminates  when  any  node  appears  in  the  labeled  set  for  both 
trees.  An  interesting  feature  of  this  procedure  is  that  the  node  which  first  appears  in 
both  trees  may  not  be  present  in  the  shortest-path,  but  enough  information  has  been 
developed  to  find  the  shortest-path.  The  speed  improvement  of  this  new  algorithm 
results  from  the  empirical  observation  that  many  fewer  nodes  need  to  be  scanned 
under  the  new  scheme.  The  new  scheme  is  more  complicated  to  implement,  but 
results  in  a  substantially  faster  software  package. 

Publication  Status 

This  paper  has  been  accepted  for  publication  in  Computational  Optimization  and 

Application- 


Title 

A  Direct  Simplex  Algorithm  for  the  Network  Flow  Problem  with  Piecewise  Linear 
Costs 

(January  1993) 

Authors 

Rajluxmi  V.  Murthy 
Richard  V.  Helgason 

Executive  Summary 

This  paper  presents  a  specialization  of  the  simplex  algorithm  for  a  network  problem 
having  a  cost  function  which  is  piecewise  linear  and  convex.  The  basis 
characterization  is  the  same  for  this  problem  but  an  improvement  can  be  made  in  the 
simplex  pricing  operation.  The  specialization  results  in  at  least  a  50%  speed 
improvement  over  existing  methods  for  this  problem. 

Publication  Status 

This  paper  has  been  submitted  for  publication  and  is  under  review. 
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APPENDIX  A 


NETWORK  FLOWS 


Richard  V.  Helgason 
Jeffery  L.  Kennington 


Department  of  Computer  Science  and  Engineering 
Southern  Methodist  University 
Dallas,  Texas  75275 

15  September  1992 


1  Introduction 

Our  aim  is  to  (1)  summarize  the  ideas  fundamental  to  efficient  implementation 
of  the  primal  network  simplex  algorithm  and  its  extensions,  and  (2)  indicate  how 
these  algorithms  may  be  effectively  used  within  other  optimization  algorithms. 

1.1  Set  Notation 

For  the  most  part  we  adopt  standard  set  notation  conventions.  Sets  will  usually 
be  denoted  by  upper  case  Roman  letters  such  as  A'.  The  empty  set  will  be 
denoted  by  $  For  a  finite  set  A'  we  let  #A'  denote  the  number  of  elements 
in  A.  We  let  in  =  and  jn  ~  {0,  ...,n}.  Given  set  A,  we  define  the 

equality  relation  on  X  to  be  A==A  =  {(x,x)  :  x£A'}.  We  will  also  use  multisets 
in  which  a  repetition  factor  is  allowed  for  set  elements.  For  a  finite  multiset  Y 
then,  #>'  will  also  incorporate  multiplicities. 

1.2  Matrix  and  Vector  Notation 

Matrices  will  usually  be  denoted  by  upper  case  Roman  letters  such  as  A.  Row 
vectors  will  usually  be  denoted  by  lower  case  Greek  letters  such  as  t.  Column 
vectors  will  usually  be  denoted  by  lower  case  Roman  letters  such  as  x.  The 
element  in  the  ith  row  and  jth  column  of  matrix  A  will  be  denoted  by  Aij. 
The  row  (column)  vector  whose  entries  are  from  the  ith  row  (jth  column)  of 
A  will  be  denoted  by  Ai  {A  j).  The  ith  element  of  a  vector  such  as  x  will 
be  denoted  by  Xj.  We  will  allow  extensive  subscripting  and  superscripting  of 
matrices  and  vectors  for  identification  purposes.  Inasmuch  as  this  may  interfere 
with  the  subscripting  convention  for  element  identification  above,  we  also  adopt 


the  functional  notation  (•)*  and  (  ),j  for  vector  and  matrix  element  identification, 
respectively,  so  that  (A\j)P?  is  the  p^th  element  of  matrix  Xt] .  We  will  use  ei 
(e*)  for  the  column  (row)  vector  whose  ith  element  is  a  1  and  whose  other 
elements  are  all  zeros.  We  will  use  e,y-  to  denote  the  column  vector  whose  ith 
element  is  a  1,  whose  jth  element  is  a  —1,  and  whose  other  elements  are  all 
zeros,  so  that  is  —  ej.  We  will  use  0  and  1  as  row  or  column  vectors  with 
orientation  and  dimension  given  by  context,  having  as  uniform  elements  0  or  1, 
respectively.  We  abuse  notation  by  allowing  in  =  {1, . . . ,  n)  to  also  be  used  as  a 
row  vector.  The  diagonal  of  matrix  A  is  the  set  of  elements  {Aj,}.  The  matrix 
A  is  said  to  be  upper  (lower)  triangular  if  A,-;  =  0  when  j  >  i  (t  >  j)  and,  more 
simply,  triangular  in  either  case.  The  matrix  A  is  said  to  be  diagonal  if  it  is  both 
upper  and  lower  triangular.  A  triangular  matrix  will  be  nonsingular  when  its 
diagonal  elements  are  all  nonzero.  The  matrix  A  is  said  to  be  triangulanzable 
if  it  can  be  brought  to  nonsingular  triangular  form  by  a  sequence  of  row  and 
column  interchanges. 

1.3  Graph  Notation 

We  define  a  set  of  nodes  or  vertices  V  to  be  any  set  of  consecutive  integers  which 
we  typically  take  to  be  in  or  jn.  Given  a  set  of  nodes  V,  we  define  an  arc  or 
edge  for  V  to  be  any  ordered  pair  (i,  j)  with  i£V,  j€V,  and  i^j.  The  arc  (t,  j) 
is  said  to  be  incident  on  (touch)  both  i  and  j ,  to  connect  i  and  j  (or  j  and  i), 
and  to  be  directed  from  i  to  j.  Formally,  a  network  or  directed  graph  is  defined 
to  be  G  —<  V,  E  >  where  V'  is  a  set  of  nodes  and  E  is  a  set  of  arcs  for  V. 
Apparently  then  E  C  (V’xV')\(V=V').  When  V  =  $  then  also  E  =  4>  and  in 
this  case  G  is  called  the  trivial  graph.  We  shall  also  allow  E  to  be  a  multiset 
when  it  is  desirable  to  have  more  than  one  arc  connect  two  nodes,  in  this  case 
one  could  more  properly  refer  to  G  as  a  multigraph.  For  #£  =  m,  we  will  find 
it  convenient  to  label  the  arcs  with  elements  from  i m. 


1.4  Visual  Representation 

The  nodes  of  a  network  may  be  viewed  as  locations  or  terminals  where  a  given 
commodity  can  be  moved  from,  to,  or  through  and  the  arcs  of  a  network  may  be 
viewed  as  unidirectional  means  of  commodity  transport  connecting  or  serving 
those  nodes.  Hence  arcs  may  represent  streets  and  highways  in  an  urban  trans¬ 
portation  network,  pipes  in  a  water  distribution  network,  or  telephone  lines  in 
a  communication  network.  The  structure  of  the  network  can  be  displayed  by 
means  of  a  labeled  drawing  in  which  nodes  are  represented  by  circles  or  boxes 
and  arcs  are  represented  by  line  segments  incident  on  two  nodes.  Each  line  seg¬ 
ment  will  have  an  arrowhead  placed  somewhere  on  it  to  indicate  the  direction  of 
the  associated  commodity  transport.  Typically  the  arrowhead  will  be  incident 
on  the  node  to  which  the  commodity  is  being  transported.  An  example  network 
illustration  is  given  in  Figure  1. 


Figure  about  here 
Figure  1:  Example  network 

1.5  Node- Arc  Matrix  Representation 

The  structure  of  a  network  may  also  be  described  using  a  node-arc  incidence 
matrix  A  given  by 


if  arc  k  is  directed  away  from  node  i, 
if  arc  k  is  directed  toward  node  i, 
otherwise. 


Apparently  then  A  *  =  e^  for  some  i  and  j,  and  we  shall  allow  ourselves  to  abuse 
notation  by  saying  that  in  this  case  the  Arth  arc  is  e,;.  An  example  node-arc 
incidence  matrix  corresponding  to  Figure  1  is  given  below. 

arcs 


1 

nodes  2 

3 

4 


1  2  3  4  5 

110  0  0 
l  -1  1  1  0 

00-101 
0  0  0  -1  -1 


6 

0 

0 

-1 

1 


7 

-1 

0 

0 

1 


1.6  Subgraphs 


A  graph  G'  =<  V',E’  >  is  said  to  be  a  subgraph  of  G  ~<  V,  E  >  if  V  C  V 
and  E'  C  E .  Note  that  G'  is  required  to  be  a  graph  itself,  so  that  V'  and  E' 
%  cannot  simply  be  arbitrary  subsets  of  V  and  E,  respectively  Further,  G'  is 

said  to  span  G  or  G'  is  said  to  be  a  spanning  subgraph  for  G  when  E'  =  E. 
Given  a  node  subset  V'  C  V,  we  define  the  subgraph  gentratiJ  by  V'  to  be 
G(V')  —  {(«,i)  6  G  :  i  G  V'  and  j  G  V').  Example  subgraphs  corresponding 
to  Figure  1  are  given  in  Figure  2. 


Figure  about  here 


Figure  2:  Example  network  subgraphs 


1.7  Paths  and  Cycles 

Given  a  graph  G  —<  V,  E  >,  a  finite  odd  length  sequence 


P-  {^i.  e«iji  i  •  •  •  i  Uj ,  e,lJt ,  } 
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whose  odd  elements  are  nodes  of  V  and  whose  even  elements  are  arcs  of  E 
is  defined  to  be  a  walk  in  C  of  length  q  >  0  in  G  if:  (1)  P  has  at  least  one 
node,  and  (2)  for  0  <  r  <  q,  arc  el(.Jr  connects  ty  and  ty+1.  Apparently  then 
from  (2),  titjT  could  be  either  (ty,  vr+i)  or  (ty+i,ty).  The  sequence  formed  by 
reversing  the  elements  of  P  is  also  a  walk  and  will  be  denoted  by  rev{P).  If  we 
envision  moving  from  v\  to  e,+  i,  utilizing  the  sequential  elements  of  P  in  order, 
we  can  assign  an  (implied)  orientation  to  the  arcs  in  the  walk  by  defining  the 
onentation  function 


0(,.  .  \  -  J  +  1  if  ei'E  =  K,fr+l)> 
iteirj,=(Vr  +  uVr). 


If  the  sequence  of  nodes  {uj , . . . ,  t??+i }  from  P  is  composed  of  distinct  nodes, 
the  walk  P  is  said  to  be  a  ( simple )  path  which  links  t>i  to  t>,+j.  It  follows  that 
the  arcs  of  a  path  are  distinct.  Apparently  then  rev(P)  is  also  a  path  which 
links  to  Vi .  It  also  follows  that  any  walk  P  of  length  0  is  a  path  which  links 
vi  to  itself.  If  the  walk  P  (1)  is  of  length  at  least  two,  (2)  {ui,e,,y, , . . . ,  tt}  is 
a  path,  (3)  {t>2,  e,3j3, . . . ,  is  a  path,  and  (4)  t’i  =  t/,+1 ,  the  walk  P  is  said 
to  be  a  cycle.  Example  walks  in  the  graph  of  Figure  1  are  given  in  Figure  3. 

Figure  about  here 

Figure  3:  Walks  in  the  example  network 

Given  a  cj  le  P  it  is  possible  to  form  other  cycles  using  the  sequential 
elements  of  P  in  wrap-around  order,  i.e.,  starting  at  vm  we  can  define  a  cycle 

P  —  {•'mi  eim}m  i  vm  +  l  t  •  •  •  i  vf  >  1  eiiii '  v2i  •  •  •  .  } 

which  retains  the  essential  arc  and  node  orders  and  arc  orientations  of  P  when 
we  envision  moving  from  vm  to  vm  on  P' ■  Thus  we  will  consider  cycles  such  as 
P  and  P'  to  be  the  same  cycle  and  also  refer  to  any  of  this  set  of  equivalent  rep¬ 
resentations  as  u  cycle  on  nodes  •  The  arcs  of  a  cycle  are  generally 

distinct,  except  for  two  special  cases  which  can  arise  when  considering  cycles  of 
length  two.  Cycles  of  the  form 

{wi ,  (vi,  V2),  w?,  (wi,  W2),  Vi} 


{t'i,(^2,vi)1u2,(t;2,ui),t;i) 

do  have  ar's  which  are  not  distinct  and  will  be  called  inadmissable  cycles.  All 
other  cycles  (which  have  distinct  arcs)  will  be  called  admissable.  Apparently 
then  if  P  is  an  admissable  cycle  on  nodes  (vi, . . . ,  then  rev(P)  is  a  distinct 
cycle  on  nodes  {«!,...,  uv}  .  A  graph  G  in  which  no  admissable  cycle  can  be 
formed  is  said  to  be  acyclic  and  is  also  said  to  contain  no  cycles. 


» 


» 


» 


» 


1.8  Connectedness  and  Components 

A  graph  G  —<  V,E  >  is  said  to  be  connected  if  any  two  vertices  u  and  v  can 
be  linked  by  a  path  in  G.  The  maximal  connected  subgraphs  of  G  serve  to 
partition  G  and  are  called  components  of  G.  If  G'  =<  (t>},4>  >  is  a  component 
of  G,  v  is  said  to  be  an  isolated  node  of  G. 

1.9  Trees 

A  nontrivial  connected  acyclic  graph  is  called  a  tree.  A  graph  which  consists  of 
an  isolated  node  only  will  be  called  a  trivial  tree.  A  graph  whose  components 
are  all  trees  is  called  a  forest.  An  endnode  of  a  tree  is  a  node  which  has  only 
one  arc  of  the  tree  incident  on  it.  A  leaf  of  a  tree  is  an  arc  of  the  tree  incident 
on  an  endnode. 

A  tree  G  —<  V,  E  >  has  several  important  properties: 

(1)  E  has  one  less  arc  than  V'  has  nodes,  i.e.  #E  =  #V'  -  1, 

(2)  if  an  endnode  and  a  leaf  incident  on  it  are  removed  from  G,  the 
resulting  subgraph  is  also  a  tree, 

(3)  if  G  has  an  arc  ( i,j )  incident  on  two  endnodes,  then  V  —  {i,  j}  and 
£  =  {(,,»}  or  £={(;,»)}. 

(4)  if#£  =  l,Ghas  exactly  one  leaf, 

(5)  if  #£  >  1,  G  has  at  least  two  leaves, 

(6)  for  every  distinct  pair  of  nodes  u,v€£,  there  is  a  unique  path  in  G 
linking  ti  to  v. 

An  example  tree  is  given  in  Figure  4. 


j  Figure  about  here 
Figure  4:  Example  tree 

A  root  is  a  node  of  a  tree  which  we  wish  to  distinguish  from  the  others, 
usually  for  some  algorithmic  purpose.  Occasionally  this  may  be  made  explicit 
by  drawing  a  tree  oriented  with  the  root  at  the  top  of  the  diagram.  Alternatively, 
this  may  be  made  explicit  in  a  diagram  by  drawing  a  short  line  segment  (with 
or  without  an  arrowhead)  incident  on  only  the  root  node.  An  example  rooted 
tree  is  given  in  Figure  5. 


A- 5 


Figure  about  here 


Figure  5:  Example  rooted  tree 


1.10  Tree  Solution  Algebra 

Consider  the  solution  of  the  system 


Ax  —  b. 


(1) 


where  A  is  the  nx(n  —  1)  node-arc  incidence  matrix  for  a  tree  T  =<  V,  E  > 
with  n  nodes  and  n  -  1  arcs  and  6  is  a  given  n-vector.  A  procedure  which  can 
be  used  to  reduce  the  system  to  a  permuted  diagonal  form  is  given  below: 


procedure  DIAGONAL  REDUCTION 


inputs:  T  =<  V',  E  > 

P 

A 

b 


nontrivial  tree  with  n  nodes 
root  node  for  T 

node- arc  incidence  matrix  for  T 
n-vector 


output:  Ax  -b  -  diagonalized  system  equivalent  to  Ax  —  b 

begin 

[Initialize] 

V  <=  U,  E  <=  E,  A  •<=  A,  6  <=  6  ; 

[Iterate  until  the  tree  becomes  trivial] 
while  E  do 
[Pick  a  leaf] 

select  an  endnode  p  not  the  root  p  in  the  tree  T  =<  V,E  >. 

let  r  be  the  other  node  of  the  leaf  incident  on  p. 

let  (i,  j)  be  the  leaf  incident  on  p  and  r. 

let  c  be  the  column  of  A  corresponding  to  the  leaf  (i,j)  ■ 

[Pivot  the  system  on  the  selected  endnode] 

br  <=  br  +  bp.  A,  <=  Ar  -f  Ap  ; 
bp  <=  ApC(bp),  Ap  <=  Ape(Ap  ); 

[Update  the  tree  by  removing  the  leaf] 

endwhile 

end 


Note  that  at  each  pivot  step  in  the  above  procedure,  a  partial  sum  of  compo¬ 
nents  of  6  is  produced  in  the  component  br,  and  in  the  last  pivot,  bp  =  6* 

is  nroduced,  Also,  after  each  pivot  and  tree  update,  the  subset  of  the  rows  and 


columns  of  A  corresponding  to  the  nodes  of  V  and  the  arcs  of  E,  respectively, 
is  the  node-arc  incidence  matrix  of  the  updated  tree  T. 

In  a  node-arc  incidence  matrix  A  for  a  tree,  let  c  be  the  column  corresponding 
to  a  leaf  with  an  endnode  p  not  the  root  node  p  and  let  r  be  the  other  leaf  node. 
Row  p  contains  only  one  nonzero  in  column  c,  and  row  r  contains  the  only  other 
nonzero  in  column  c.  Thus  when  a  pivot  is  made  at  the  matrix  position  Apc, 
Arc  will  be  zeroed  and  Apc  will  become  1  if  it  was  -1,  but  no  other  entries  in 
A  will  be  altered. 

Now  A  initially  has  2n  —  2  nonzeros  and  n  -  1  pivots  are  made  overall,  so 
that  the  final  pivot  produces  a  A  with  n  —  1  nonzeros,  all  in  the  pivot  rows. 
Thus  row  p  of  A,  the  only  row  in  which  no  pivot  occurred,  must  contain  all 
zeros.  It  follows  that  the  system  {  1)  has  a  solution  if  and  only  if  bp  —  0,  hence 
no  solution  is  possible  unless  5*  =  0  .  Furthermore,  since  n  -  1  pivots 
have  been  made,  the  matrix  A  has  rank  n  -  1  so  that  when  =  0,  the 

solution  produced  by  the  algorithm  is  the  unique  solution  to  system  (1). 

To  illustrate  the  use  of  the  algorithm  consider  the  tree  in  Figure  4.  The 
original  system  corresponding  to  (  1)  is 


-1 

0 

0 

1 


0  0 

1  1 

-1  0 

0  -1 


" 

r  ^ 

'  *1  ' 

*1 

62 

*2 

— 

63 

*3 

L 

°4 

Selecting  node  4  as  the  root  and  using  the  sequence  of  selected  endnodes  3,2, 1 
produces  the  following  sequence  of  systems. 


-10  0 
0  0  1 

0  1  0 

1  0  -1 


-  . 

'  b ,  ’ 

*1 

62  +  b  3 

*2 

-63 

.  *3  . 

b 4 

O 

O 

1 

L— 

r  ^ 

0  0  1 

*i 

62  +  bz 

0  1  0 

-63 

I  0  0 

J 

*3 

62  +  63  +  64 

’10  0' 

_  _ 

0  0  1 

*1 

0  1  0 

*2 

55 

0  0  0 

„  *3  . 

-6, 

62  +  63 

-63 

6l  +  62  +  ^3  +  64 


Now  let  us  consider  adjoining  an  additional  column  e*,  where  1  <  k  <  n, 
to  the  right  of  A  and  lengthening  x  to  accomodate  the  additional  column.  The 
expanded  system  is  then 


[<4|e*]x  =  6, 


(2) 


Suppose  that  we  agree  to  choose  k  as  the  root  node  (p  =  k)  and  apply  the  above 
procedure  to  the  expanded  matrix  and  original  tree,  i.e.  in  the  initialization  step 
we  set  A  <=  [A\ep]  instead  of  A  <=  A. 

The  same  vector  6  and  matrix  A  is  produced  in  the  first  n  —  1  (original) 
columns  and  A  „  =  e*  =  ep.  The  system  (  2)  is  also  permuted  diagonal  but  is 
now  of  rank  n  and  in  its  solution  has  e*  =  2Ls=i  6*  ■  Since  [ylfe*]  is  square,  the 
solution  produced  by  the  procedure  must  be  unique.  Furthermore,  the  solution 
to  (  1)  produced  by  the  above  procedure  when  ^t=i  6*  =  0  must  have  the  same 
first  n  —  1  variable  values  as  those  produced  for  the  enlarged  system  (  2). 

In  the  above  example,  the  original  system  corresponding  to  (  2)  with  root 
node  4  is 


'  -1  0  0  0' 

p  - 

*1 

'  bi  ' 

0  110 

*2 

bi 

0-1  0  0 

*3 

63 

1  0-11 

.  z«  . 

.  . 

And  after  the  same  sequence  of  pivots,  the  final  equivalent  system  produced  is 


'  1  0  0  0  ‘ 

’  *1  " 

0  0  10 

*2 

h  +  b  3 

0  10  0 

*3 

-63 

0  0  0  1 

.  *4  . 

61+62  +  63  +  64 

We  remark  that  this  usage  of  an  extra  e *  column  in  conjunction  with  the 
solution  of  system  (  1)  provides  strong  impetus  to  extend  the  node-arc  matrix 
representation  to  include  a  representation  of  a  root  it  by  a  column  e*,  when 
the  underlying  graph  is  a  tree.  We  will  find  it  useful  to  do  so  even  when  the 
underlying  graph  is  a  tree  with  additional  arcs. 

2  Primal  Simplex  Algorithm 

All  the  network  models  presented  in  this  Chapter  may  be  viewed  as  special 
cases  of  the  linear  program  and  many  of  the  algorithms  to  be  presented  are 
specializations  of  the  primal  simplex  algorithm.  Familiarity  with  the  simplex 
method  is  assumed  and  this  section  simply  presents  a  concise  statement  of  the 
primal  simplex  method  which  will  be  specialized  for  several  of  the  network 
models. 

Let  A  be  an  mxn  matrix  with  full  row  rank,  c  and  u  be  n-component 
vectors,  and  6  be  an  m-component  vector.  Let  X  =  {x  :  Ax  =  b,  6  <  x  <  u} 
and  assume  that  A'  ^  $.  The  linear  program  is  to  find  an  n-component  vector 
x  such  thai  cx  =  min{cx  :  x^X}.  We  adopt  the  convention  that  the  symbol 
x  refers  to  the  vector  of  unknown  decision  variables  and  the  symbol  x  refers  to 
some  known  point  in  X. 

Since  A  has  full  row  rank,  by  a  sequence  of  column  interchanges  it  may 
be  displayed  in  the  partitioned  form  [5|Ar],  where  B  is  nonsingular.  Writh  a 


corresponding  partitioning  of  both  x  and  u,  we  may  write 

A  =  {(xB|xN)  :  Bxb  +  Nxn  =6,0<xb  <ub,0<xs  <un) 

A  point  (xB  I  xN)  €  A-  is  called  a  basic  feasible  solution  if  x*  €  {0,  for 
all  j  €  *(n  — rn)  .  The  variables  in  xB  ( xN )  are  called  the  basic  ( nonbastc ) 
variables.  It  is  well  known  from  linear  programming  theory  that  for  every  linear 
program  with  X  £  there  exists  at  least  one  basic  feasible  solution  which  is 
optimal. 

We  say  that  two  basic  feasible  solutions  are  adjacent  if  xBl  and  xB>  differ 
by  exactly  one  variable,  so  that  B\  may  differ  from  B2  in  exactly  one  column. 
The  primal  simplex  algorithm  is  an  iterative  procedure  which  begins  with  some 
basic  feasible  solution  and  moves  through  a  sequence  of  adjacent  basic  feasible 
solutions  until  a  stopping  criterion  is  met  which  guarantees  that  the  final  basic 
feasible  solution  is  optimal. 

A  mathematical  description  of  the  primal  simplex  algorithm  follows: 

procedure  PRIMAL  SIMPLEX 

inputs:  A  -  mxn  constraint  matrix 

c  -  m-component  vector  of  unit  costs 
u  -  n-component  vector  of  upper  bounds 
6  -  m-component  right-hand-side  vector 

output:  (xB  |  x^),  an  optimal  basic  feasible  solution  for  min{cx  :  x€A'} 


assumptions:  1.  rank(A)  =  m 

^  2.  A  =  {x  :  Ax  =  6 , 6  <  x  <  u}  ^  $ 

begin 

[Initialization] 
optimal  e=  ‘no’; 

let  (xB  lxy)  €  X  be  any  basic  feasible  solution  with  A ,  c,  and  u  parti- 
0  tioned  in  corresponding  fashion  as  [B  |  A],  (c®  | c7^),  and  (uB  |  u^); 

while  optimal  =  ‘no’  do 
[Dual  Calculation] 


*<=cbB~' 

[Pricing] 

L<=  {j  :  rN.j  -  c?  >  0  and  x^  =  0} ; 
M  <=  {j  :  xN  j  —  <  0  and  xj^  =  Uj1}  -, 
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then 

optima]  <=  ‘yes’ 

else 

select  k  6  L  U  M  ; 
[Column  Update] 


[Ratio  Test] 


y  <=  B~1Nk\ 


R  <=  {«  :  j h  >  0}  and  S  <=  {»  :  y,  <  0}; 

if  it  G  L 
then 

A  <=  min{uj  ,minie*  .min^s 
T  <=  {«'  G  R  :  (£)  =  A}  U  {i  6  S  :  (^- ) 

else 


A  <=  minfii^ 
T<={i£R  : 


(^)=A,U  VtS:(d) 


mml€fl 
..a  #b 


endif 

[Value  Update] 


ifibei 

then 


<=  A  and  xB  <=  xB  -  Ay; 

else 

x*1'  <=  -  A  and  xB  <=  xB  Ay; 

endif 
if  A/u" 
then 

[Basis  Exchange  Update] 


select  j  e  T  ; 

Interchange  xB  with  x£  to  form  a  new 
partitioning  [B  |  N]  ; 

endif 


> 


» 

endif 
end  while 
end 

)  The  only  sticky  issue  concerning  the  above  algorithm  is  that  it  may  be  pos¬ 

sible  to  move  through  a  sequence  of  basis  interchanges,  all  of  which  correspond 
to  the  same  actual  point,  only  changing  the  representations  to  correspond  to 
the  varying  bases.  That  is,  there  may  exist  a  sequence  of  basis  interchanges 
all  having  A  =  0  which  could  result  in  the  above  procedure  being  stuck  in  a 
nonterminating  loop.  This  phenomenon  is  known  as  cycling.  Much  discussion  of 
£  cycling  can  be  found  in  the  literature  along  with  anticycling  rules5  which  have 

been  developed  to  insure  convergence. 


3  Linear  Network  Models 

Let  <  V,  E  >  be  a  network  through  which  some  commodity  will  be  flowing. 
Associated  with  each  node  v  6  V  we  define  a  requirement  r„  A  node  having  a 
supply  of  the  commodity  is  assigned  a  positive  requirement  equal  to  the  supply. 
Conversely,  negative  requirements  correspond  to  demands  for  the  commodity  at 
the  specified  nodes.  The  requirements  for  all  transshipment  nodes  are  zero. 

Suppose  that  for  the  example  network  illustrated  in  Figure  1,  that  nodes  1 
and  3  are  supply  nodes  with  supply  of  10  and  5,  resp<  tively;  and  that  node  4 
is  a  demand  point  with  a  demand  of  15.  This  netwoi  and  the  corresponding 
requirements  are  illustrated  in  Figure  6. 


Figure  about  here 


Figure  6:  Example  network  with  requirements 


For  linear  network  problems,  one  is  seeking  a  set  of  flows  on  the  arcs  which 
satisfy  the  supply  and  demand  restrictions  at  each  node.  For  the  example 
network,  the  total  flow  into  node  2  plus  the  five  supply  units  which  originate  at 
node  2  must  equal  the  total  flow  departing  node  2.  We  say  that  a  feasible  flow 
satisfies  flow  conservation  which  implies  that  it  is  an  element  of  {x  :  Ax  =  r}, 
where  A  is  the  corresponding  node-arc  incidence  matrix  for  <  V,  E  >. 

1  We  have  never  observed  cycling  in  any  of  our  software  implementations  of  this  algorithm 
and  have  not  incorporated  any  of  the  anticycling  rules  in  our  software. 
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For  the  example  network,  the  flow  conservation  equations  are 


’  *1  ' 

■  1  1  0  0  0  0  -1  ‘ 

*2 

'  10  ' 

-1-111000 

*3 

5 

0  0-1  0  1-1  0 

*4 

= 

0 

0  0  0  -1  -1  1  1 

*5 

-15 

* 

*6 

.  17  . 

Associated  with  each  of  the  arcs  in  E,  we  define  a  vector  of  unit  flow  costs  c 
and  a  vector  of  flow  capacities  or  bounds  u.  Thus,  the  cost  for  each  unit  of  flow 
in  arc  k  is  given  by  c*  and  the  flow  on  arc  k  is  restricted  to  the  interval  [0,u*]. 
Mathematically,  the  linear  network  model  on  <  V,  E  >  with  node-arc  incidence 
matrix  A  is  given  by 

min{cx  :  Ax  =  r,0  <  x  <  u} 


A  sample  model  for  the  network  illustrated  in  Figure  6  is 


minimize  xj  +  3x2  +  5x3  -  7x4  +  7xj  —  X6  +  9xr 
subject  to 
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3.1  Basis  Characterization 

Recall  that  one  of  the  assumptions  for  the  primal  simplex  algorithm  is  that 
the  constraint  matrix  has  full  row  rank,  i.e.  the  mxn  constraint  matrix  A 
has  rank(A)  =  m.  But  it  is  clear  that  if  A  is  a  node-arc  incidence  matrix, 
rank(A)  <  m  since  for  a  connected  graph,  the  rank  of  a  node-arc  incidence 
matrix  is  one  less  than  the  number  of  nodes  (also  \A  —  0).  Traditionally,  a  root 
arc  is  appended  to  the  problem  so  that  the  constraint  matrix  will  have  full  row 
rank.  Let 
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A'£  =  {(x|a)  :  (/l|e,] 


a 


The  revised  model  is  then  simply 

min{cx  :  (x|a)  6  XE], 

where  the  enlarged  constraint  matrix  [A|e,]  has  full  row  rank,  so  that  the  primal 
simplex  algorithm  can  be  applied  directly  to  this  model.  In  Section  2.10  it  was 
shown  that  the  root  arc  will  carry  no  flow  (a  =  0)  when  JV  r*  =  0  and  we  are 
only  solving  the  linear  system.  This  will  also  be  true  when  any  nonbasic  arcs 
are  set  to  upper  bound,  since  setting  a  nonbasic  flow  xk  =  ut,  where  x*  is  the 
flow  on  arc  (i,j)  is  equivalent  to  adding  uk  to  r j  and  subtracting  uk  from  r,, 
thus  preserving  the  condition  r,  =  0. 

Let  B  be  a  basis  for  [A\ep]  so  that  the  entire  matrix  is  partitioned  as  [J5| A'']. 
B  may  be  further  partitioned  into  [S|e,|.  Recall  that  the  arcs  corresponding  to 
S  must  form  a  spanning  tree  from  <  V,E>.  We  may  write  the  corresponding 
basic  solution  as  (xs|aJxN). 

By  row  and  column  interchanges  B  may  be  displayed  in  lower  triangular 
form.  The  trees  corresponding  to  the  bases 
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are  illustrated  in  Figure  7. 


Figure  about  here 

Figure  7:  Basis  trees  for  sample  network 

The  trees  can  be  used  to  determine  the  row  and  column  interchanges  required 
to  display  the  bases  in  lower  triangular  form.  The  root  node  and  root  arc  are 
always  placed  in  the  last  row  and  column.  The  first  m—  1  rows  and  columns  are 
determined  recursively  by  a  process  in  which  we  first  select  an  endnode  with  its 
corresponding  leaf  arc  and  then  remove  them  from  the  tree  after  appropriately 
reordering  the  matrix.  Note  that  this  ordering  may  not  be  unique  since  multiple 
endnodes  are  always  present, 

A  procedure  which  can  be  used  to  display  a  linear  network  basis  in  lower 
triangular  form  is  given  below: 


I 
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procedure  DISPLAY  LOWER  TRIANGULAR 


input: 


-  mxm  basis  for  a  linear  network  problem 


outputs:  norder\i ]  -  the  node  in  row  position  i 

arcorder\j]  -  the  arc  in  column  position  j 

assumption:  S  corresponds  to  a  spanning  tree 

begin 

[Initialization] 

let  <  V,  T  >  be  a  network  corresponding  to  5 
nodeorder[m ]  <=  p,  arcorder[m]  <=  ep; 

» <=  1; 

while  i  <  m  do 
[TVee  Reduction] 

let  v  6  V  be  an  endnode  of  <  V,T  >; 
let  e;i  be  the  leaf  arc  incident  to  v; 

nodeorder[i]  e=  v ,  arcorder\i]  <=  e^,  V  <=  T  <=  T\{e jk},  i  <=  «  +  l; 

endwhile 

end 


An  app'  cation  of  procedure  DISPLAY  LOWER  TRIANGULAR  to  the  tree 
in  Figure  7a  with  certain  choices  in  the  tree  reduction  yields  nodeorder  = 
[3, 1,4,2]  and  arcorder  —  [e^s,  e4i, €24,62).  The  corresponding  matrix  is: 


(nodes) 

3 
1 

4 

2 


e2<  *2 

0  0 

0  0 

-1  0 

1  1 


Another  application  of  the  procedure  with  different  choices  made  in  the  tree 
reduction  would  yield  nodeorder  =  [1,3, 4, 2)  and  arcorder  =  [e*!, e23,e24, ej], 
with  corresponding  matrix  is: 


(nodes) 

1 

3 

4 
2 


«24 

0  0 
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3.2  Dual  Calculation 


Based  on  the  results  in  Section  3.2,  it  is  known  that  every  basis  for  a  linear 
network  problem  takes  the  form  [Sje^j,  where  the  columns  of  S  correspond  to  a 
spanning  tree  <  V,  T  >  and  S  is  triangularizable.  Therefore  the  dual  calculation 

r[S\efi)  =  [cs\0] 

can  be  specialized  to  exploit  the  underlying  network  structure.  Since  S  consists 
of  columns  from  a  node-arc  incidence  matrix,  every  column  of  5  is  a  vector  e^ 
for  some  *  and  j.  Hence  5  =  [e,Ul,  e,^,, . . .  ,eim_,jm_,}  and  the  dual  calculation 
reduces  to  the  system 


Since  this  system  is  upper  triangular  it  can  be  solved  by  back  substitution. 
Beginning  with  rrp  =  0,  the  duals  for  all  basic  arcs  incident  on  node  p  can  be 
determined  Once  these  duals  are  known,  all  duals  for  basic  arcs  incident  on 
those  arcs  can  be  determined.  Continuing  in  this  manner,  eventually  all  duals 
will  be  determined.  The  following  procedure  formally  states  how  this  can  be 
accomplished  without  actua  ly  triangularizing  a  matrix. 

procedure  DUAL  CALCULATION 

inputs:  [S|ep]  -  mxm  basis  for  a  linear  network  problem 

cs  m  —  1  vector  of  basic  costs,  where 

cfj  is  the  unit  cost  for  basic  arc  el; 


ou'  put:  x  -  *•  =  [cs|0][S|e^]_1 


begin 

[Initialization] 

let  <  V,  T  >  be  the  spanning  tree  corresponding  to  5; 
for  i  =  1, . . m 
/a6e/[i]  =  ‘no’; 

endfor 

*p  <=  0,  labe![p]  <=  ‘yes’,  k  <=  1; 


while  k  <  m  do 

let  tij  £  T  such  that  /a6e/[i]  /  label[j}\ 
if  lat>el[j]  =  ‘yes’ 

then 

-  cfj,label[j]  =  'yes’; 

else 

*i  =  *j  ■ f  c5,/a6e/[t]  =  ‘yes’; 

endif 

k  <=  k  +  1; 
endwhile 
end 


Consider  the  basis  tree  in  Figure  8  with 

lc12>  cf 4,  c4l>  C51  >  c74  ’  c76]  =  (1, 2,3,  —  1,  —  2,  —3] 
Since  node  4  is  the  root  node,  —  0.  The  equations 


*3  -  »4  =  2  1 

*4  —  *i  =  3  > 

*7  ~  *4  =  ~2  J 


yield  x-3  =  2,  jtj  =  -3,  and  ■K^  =  -2.  Then  the  equations 


*5  “  *1  =  -  1  1 

*1  ~  *2  =  1  > 

*7  “  *6  ~  —3  J 


yield  rr5  =  -4,jt2  =  -4,  and  T6  =  1.  Also  note  that  the  dual  calculation  only 
involves  addition  and  subtraction.  Hence,  if  the  cost  coefficients  are  all  integer, 
then  the  dual  variables  will  also  be  integral  valued. 


Figure  about  here 


Figure  8:  Dual  variable  calculation  example 


3.3  Column  Update 

Suppose  e„  denotes  the  nonbasic  arc  which  prices  favorably  and  is  selected  for  a 
potential  flow  change.  There  is  no  guarantee  that  the  flow  will  actually  change 
since  the  ratio  test  could  yield  A  =  0.  The  column  update  step  of  the  primal 
simplex  algorithm  requires  that  y  =  [5je,]-1e,r  be  determined.  Since  the  arcs 
of  5  correspond  to  a  spanning  tree  and  [Sje,,]  is  lower  triangular,  the  calculation 
of  y  can  be  simplified. 


The  updated  column  y  is  a  vector  which  solves  the  lower  triangular  system 
[5|ep]y  =  e,t  or,  in  component  form, 

S.iyi  +  S  2y2  4-  •  ■  •  4-  S  m-l  ym  —  1  +  tpVm  — 

That  is,  a  set  of  scalars  jn , . . . ,  j/m  are  sought  which  when  multiplied  by  the 
columns  of  5  and  the  vector  ep  and  added  together  form  the  vector  e,« .  Let 
<  V,T  >  denote  the  tree  corresponding  to  S  and  let 

P  —  » e* j , . . . , Vq , 

denote  the  simple  path  in  <  V,T  >  linking  t  to  t.  Such  a  path  is  illustrated 
in  Figure  9.  Since  some  arcs  in  the  path  may  be  directed  uum  some  t>,  toward 
and  others  from  some  toward  t-;,  arrow  heads  have  been  omitted  from 
the  illustration. 


Figure  about  here 

Figure  9:  A  simple  path  linking  «  to  t 


By  reordering  the  rows  and  the  arcs  in  5,  the  system  of  equations  corre¬ 
sponding  to  the  arcs  in  the  path  may  be  illustrated  as  shown  below. 
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Since  the  first  row  of  the  above  system  has  a  single  nonzero,  yi  is  uniquely 
determined  and  must  be  either  1  or  —1.  Once  yi  is  known,  j/2  is  then  uniquely 
determined  and  must  be  either  1  or  -1.  Similarly,  y3, . . . ,  yf  can  be  determined 
successively.  Therefore,  a  solution  to  [S|ep]y  =  e*«  can  be  constructed  by  setting 
the  components  of  y  corresponding  to  the  path  from  s  to  t  in  <  V,  T  >  to  ±ls 
as  described  above  and  setting  all  other  components  to  zero  This  is  generally 
called  a  cycle  trace  and  is  formally  presented  in  the  following  procedure. 

procedure  CYCLE  TRACE 

inputs:  [5|ep]  -  mx m  basis  for  a  linear  network  problem 

ett  -  m  vector  corresponding  to  a  selected 

nonbasic  arc 


output: 


V 


V=  [S|e,]-1e.« 


begin 

[Initialization] 

let  <  V,T  >  be  the  spanning  tree  corresponding  to  5; 
let  P  = 

be  the  simple  path  in  <  V,  T  >  linking  node  s  to  node  <; 
for  i  =  l,...,Tn 
!/.  =  0; 
endfor 
k  <=  1; 

while  k  <  q  do 

let  c  be  the  column  index  of  S  corresponding  to  arc  euyt; 

e*kjk  ~  evu  ~ 

then 

Vc  =  1; 

else 

yc  =  -1; 
endif 
k  <=  k  +  1; 
endwhile 
end 

Consider  vhe  basis  tree  in  Figure  8  and  suppose  r  =  5  and  1  =  6.  The 
simple  path  and  corresponding  values  of  the  updated  column  are  illustrated  in 
Figur  10.  Note  that  the  nonzero  components  of  the  updated  column  y  are 
identical  to  the  orientation  function  (see  Section  1.7)  on  the  simple  path  from 
s  to  t.  Furthermore,  the  components  of  y  are  from  {1,-1}. 


[  Figure  about  here 
Figure  10:  Updated  column  example 


3.4  Basis  Exchange  Update 

As  seen  in  Sections  3.2  and  3.3,  the  key  operations  for  the  primal  simplex 
algorithm  can  be  performed  directly  on  the  spanning  tree  <  V,T  >.  In  this 
section  a  data  structure  used  to  store  the  spanning  tree  in  computer  memory  is 
presented  along  with  an  algorithm  vhich  will  perform  the  basis  exchange  update 
using  this  data  structure. 

Suppose  the  rooted  tree  is  drawn  in  the  plane  placing  the  root  node  at  the  top 
with  the  branches  extending  downward  as  illustrated  in  Figure  11a.  One  may 


imagine  tracing  a  line  around  the  contours  of  the  tree  as  illustrated  in  Figure  11b. 
Traversing  a  tree  in  this  way  has  become  known  as  a  depth-first  search .  For  the 
example,  the  nodes  in  this  search  could  be  ordered  as  4, 3,4, 1,5, 1,2, 1,4, 7, 6, 7, 4.  By 
eliminating  all  duplicate  occurances  an  ordering  known  as  preorder  is  obtained. 
For  this  example,  the  corresponding  preorder  is  4, 3, 1,5, 2, 7 ,6.  The  label  which 
gives  the  next  node  in  the  preorder  for  node  v  is  known  as  the  thread ,  denoted 
by  The  thread  for  the  example  is  illustrated  in  Figure  11c. 


Figure  about  here 

Figure  11:  Illustration  of  thread  labels  for  sample  rooted  tree 

The  other  two  labels  are  related  to  the  path  from  a  given  node  v  to  the  root 
p  in  the  basis  tree  <  V,  T  > .  Let 


P  —  {^1 1  eiiji  1  e»3j3i  •  ■  •  >  v«+l} 

denote  the  simple  path  in  <  V,  T  >  which  links  v  to  p.  The  predecessor  label 
p(v)  is  Vi  (the  second  node  on  the  path)  and  the  distance  label  d(v )  is  q  (the 
number  of  arcs  on  the  path).  Both  the  predecessor  and  the  distance  labels  of  the 
root  p  are  defined  to  be  zero.  The  labels  for  the  rooted  tree  shown  in  Figure  1 1 
are  given  in  Table  1. 

Table  1.  Labels  for  the  rooted  tree  in  Figure  11a 


node 

V 

predecessor 

P(v) 

thread 

<(») 

distance 

d(v) 

1 

4 

5 

1 

2 

1 

7 

2 

3 

1 

1 

1 

4 

0 

3 

0 

5 

1 

2 

2 
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4 

2 

7 
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1 

Figure  about  here 

Figure  12:  New  rooted  tree  after  basis  exchange  update 

Both  the  dual  calculation  (Section  3.2)  and  the  column  update  (Section  3.3) 
can  be  easily  implemented  in  software  using  the  triple  labels  to  represent  the 
basis  tree  <  V,T  >.  The  only  tricky  part  is  the  technique  needed  for  a  basis 
exchange  update.  For  example,  suppose  the  arc  (2,7)  is  exchanged  with  the  arc 
(1,4)  in  the  rooted  tree  illustrated  in  Figure  11. 


Table  2.  Labels  after  the  basis  exchange  update 


node 

V 

predecessor 

pH 

thread 

f(r) 

distance 

d(v) 

1 

2 

5 

3 

2 

7 

1 

2 

3 

1 

7 

1 

4 

0 

3 

0 

5 

1 

6 

4 

6 

7 

4 

2 

7 

4 

2 

1 

The  new  tree  is  illustrated  in  Figure  12  and  the  corresponding  triple  labels 
are  given  in  Table  2.  An  algorithm  which  will  perform  this  exchange  is  given 
below. 

procedure  BASIS  EXCHANGE  UPDATE 


inputs: 

<V.T> 

the  current  basis  tree 

pW 

the  predecessor  label  oi  node  i 

m 

the  thread  label  of  node  i 

d[i] 

the  distance  label  of  node  i 

(«.») 

the  arc  selected  to  become  part  of  the  new  basis 

k 

node  such  that  k  and  p[it]  are  the  end  nodes  of  the 
arc  which  is  to  be  removed  from  the  current  basis 

outputs: 

P[»] 

updated  predecessor  label  of  node  i 

no 

Updated  thread  label  of  node  i 

d[t] 

updated  distance  label  of  node  i 

assumption:  k  and  p[£]  are  on  the  path  from  u  to  the  root  node  of  the 

current  basis  tree 


begin 

[Initialization] 

q,  j  <=  p[g],  k'  <=  p[Jt],  /  <=  d[q']  +  1,  end  <=  ‘no’; 
while  end  =  ‘no’  do 

I'  <=  d[ {],  m  <=  1  -  rf[i],  d(i]  <=  /,  z  <=  *,  z  <=  <[i]; 
while  </[r]  >  /'  do 
d(x]  d[z]  +  m,  z  <=  x,  z  ^  <[*]; 
endwhile 
r  <=  >; 

while  f[r]  ^  i  do 


r  <=  <[r]; 
endwhile 

t[r]  <=  x\ 

i£i  =  k 

then 

<=  ffal,  <[«']  <=  9.  pfe]  <=  end  <=  ‘yes'; 

else 

<[2]  <=  i.  r  <=  «-  «  ■<=  j>  i  <=  p[>],  p[*3  <=  r,  /  <=  /  +  1; 
endif 
endwhile 
end 

In  our  software  implementation  of  the  primal  simplex  algorithm  for  linear 
network  problems,  the  basis  exchange  update  is  also  integrated  with  the  value 
update.  In  addition,  since  the  set  of  nodes  whose  dual  values  change  and  the 
set  of  nodes  whose  distance  labels  change  is  the  same  set,  the  dual  update 
calculation  is  also  integrated  with  the  basis  exchange  update.  All  of  these  spe¬ 
cializations  can  result  in  software  which  can  execute  one  hundred  times  faster 
than  general  purpose  linear  programming  software.  Most  of  the  modern  com¬ 
mercial  linear  programming  systems  now  contain  a  specialized  network  solution 
module. 

4  Generalized  Networks 

Let  G  be  an  mxn  matrix  with  full  row  rank  such  that  every  column  of  G  has 
at  most  two  nonzero  entries.  Let  c  and  u  be  n-component  vectors  and  b  be 
an  m-component  vector.  Let  Y  =  {x  :  Gx  =  6,  6  <  x  <  u}  and  assume  that 
y  0.  The  generalized  network  problem  is  to  find  an  n-component  vector  x 
such  that 

cx  =  minfej* :  *  €  Y) 

Note  that  the  linear  network  model  is  a  special  case  of  the  generalized  network 
problem.  The  revised  model  developed  in  Section  3.1  is  simply 

min{cz  :  (ar,o)  G  XE } 

X 


where 

XE  =  {(*l«)  =  Ax  +  epa  -  r,0  <  x  <  ti} 

The  corresponding  generalized  model  has  G  —  [/4|e^]. 

Let  V  =  tm  —  {1,.  Then  n  arcs  on  the  node  set  V  cam  be  constructed 

using  the  n  columns  of  G.  If  G  j  has  two  nonzero  entries  in  rows  i  and  k,  we  let 
the  corresponding  arc  be  (i,  it).  If  G  j  had  only  a  single  nonzero  entry  in  row  i, 
we  let  the  corresponding  arc  be  a  root  (»). 


Consider  the  following  matrix. 


120  0000  00 

-103-2000  00 

000  0122  00 

000  0040-10 

020  1000  00 

000  0001  03 

000  0000  -1  1 

Using  the  rule  presented  above,  the  following  network  can  be  constructed: 

<  {1.2, 3, 4,5, 6, 7},  {(1,2),  (1,5), (2), (2, 5), (3),  (3, 4),  (3, 6), (4. 7), (6. 7)}  > 

A  visual  representation  of  this  network  is  given  in  Figure  13.  The  numbers 
assoc’  ’ted  with  the  ends  of  the  arcs  correspond  to  the  nonzero  entries  in  the 
columns  of  G. 

Figure  about  here 

Figure  13:  Example  generalized  network 


4.1  Basis  Characterization 

Let  Gbeanmxn  matrix  with  full  row  rank  such  that  every  column  of  G  has  at 
most  two  nonzero  entries.  Let  B  be  a  basis  for  G  and  let  B  be  a  reordering  of 
rows  and  columns  of  B  such  that  B  is  displayed  in  block  diagonal  form.  That 
is, 


For  example: 


and 


B  = 


B1 


B2 


B f 


B  = 


1  2  0  0  0  0  0 

-10-2000  0 
0  0  0  1  2  2  0 

0  0  0  0  4  0  -1 

0  2  1  0  0  0  0 

0  0  0  0  0  1  0 

00  0000-1 
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with  p  =  2.  The  number  of  blocks,  p,  can  be  as  many  as  m.  Let  the  networks 
associated  with  each  of  the  blocks  be  denoted  by  <  Vl,Ex  Vp,Ep  >■ 

For  the  above  example.  <  Vl.El  >=<  {1,2,5}, {(1,2), (2,5), (1,5)}  >  and 
<  \ /2,E2  >=<  {3, 4, 6, 7},  {(3), (3, 4), (3, 6),  (4, 7)}  >.  The  corresponding  net¬ 
works  are  illustrated  in  Figure  14 

Figure  about  here 

Figure  14:  Basis  for  a  generalized  network 

A  connected  network  having  exactly  one  cycle  is  traditionally  called  a  one- 
tree  and  a  connected  network  having  exactly  one  root  arc  is  traditionally  called 
a  rooted  tree.  It  is  well  known  that  the  connected  components  from  a  basis  of 
a  generalized  network  are  either  one-trees  or  rooted  trees.  That  is,  <  Vr  1 ,  £7*  > 

. <  VP.EP  >  are  either  one-trees  or  rooted  trees.  For  the  basis  illustrated 

in  Figure  14.  one  component  is  a  rooted  tree  and  the  other  is  a  one-tree. 

4.2  Dual  Calculation 

Since  the  basis  B  has  special  structure,  irB  =  c  can  be  specialized  to  exploit 
this  structure.  Since  B  is  block  diagonal  we  obtain 


rr1  |  x2  |  ■  ■  |  xp 


=  [  C1  |  cJ  I  •  •  ■  I  <*  } 


Hence,  the  p  systems  x*B?  =  c q  cam  be  solved  independently. 

Two  cases  must  be  considered.  Suppose  Bq  corresponds  to  a  rooted  tree. 
Then  by  row  and  column  interchanges,  Bq  can  be  displayed  in  lower  triangular 
form  and  x?B?  =  c?  can  be  solved  by  back-substitution. 

Consider  the  example  illustrated  in  Figure  15  with 
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Figure  about  here 

Figure  15:  Rooted  tree  component  of  a  basis  for  a  generalized  network 

From  the  fourth  equation,  ir3  =  -2.  Then  from  the  third  equation,  *4  =  2. 
The  first  two  equations  yield  ir6  =  3  and  ir7  =  -4.  For  the  rooted  tree  case,  a 
procedure  similar  to  the  one  given  in  Section  3.2  can  be  developed. 

Figure  about  here 

Figure  16:  One-tree  component  of  a  basis  for  a  generalized  network 
Consider  the  example  illustrated  in  Figure  16  with 
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This  system  of  equations  for  a  cycle  is  almost  lower  triangular  (only  the  last 
column  has  a  nonzero  entry  above  the  diagonal)  and  takes  the  special  form: 

ai  bn 

b}  a2 
b2  03 


it  •  S 

L  |  I  [  bn  —  1  j  On  J 

wrhich  corresponds  to  the  cycle  illustrated  in  Figure  17.  Under  the  assumption 
that  the  above  system  has  full  row  rank,  then  the  unique  solution  is  given  by: 

ci/ai  +  u’!C2/a2  +  Wiw2c3/a3  + - 1-  Ujw2-  ■  tr„-i c„fa„  /ox 

TTj  _  - - - - - 

1  -  tt'lU’2-  •  Wn 


tr,+i  = - 7 - for  t  =  1, 

bi 

where  to,  =  —b,/a,  for  t  =  1, . . . ,  n . 

Figure  about  here 

Figure  17:  Cycle  for  a  generalized  network 

For  the  example  illustrated  in  Figure  17,  wj  =  l,w2  =  1/2,  w3  = 
2,  tr2  =  -3,  and  tr3  =  0. 


-l,Ti  = 


A-24 


An  extension  of  the  procedure  DUAL  CALCULATION  given  in  Section 

3.2  can  be  developed  for  the  case  of  generalized  networks.  The  basic  idea  is 
that  the  components  corresponding  to  rooted  trees  are  triangularizable  and 
the  components  corresponding  to  one-trees  are  nearly  triangularizable.  For  the 
components  corresponding  to  one-trees,  (3)  and  (4)  are  used  for  calculations 
involving  the  cycle  and  all  other  calculations  involve  a  lower  triangular  system 
of  equations.  Hence,  once  the  duals  on  the  cycle  are  known,  all  other  duals  can 
be  determined  by  back-substitution. 

4.3  Column  Update 

Suppose  that  the  nonbasic  arc  corresponding  to  Gj  prices  favorably  and  is 
selected  for  a  potential  flow  change.  The  column  update  step  of  the  primal 
simplex  algorithm  requires  that  y  be  determined  where  y  -  B~1G  j  or  By  =  G  j . 
Since  B  is  block  diagonal,  this  may  be  written  in  the  form: 


[  Bf  \  L  If  J  L  9P  J 

Therefore,  one  must  solve  p  systems  of  the  form  Bqyq  =  gq  where  the  corre¬ 
sponding  network  <  V* ,  ?q  >  is  either  a  rooted  tree  or  a  one-tree.  Furthermore. 
gq  will  have  at  most  two  nonzero  entries.  If  g q  =  6,  then  yq  —  0. 

Suppose  gq  has  two  nonzero  entries  a  and  b  in  positions  i  and  j.  respectively. 
Then  gq  —  aa  +  bej.  Let  yq  and  yl  solve  Bqyq  =  ae*  and  Bqy\  =  6e;,  respec¬ 
tively.  Then  yq  =  yl  -by*.  Hence,  a  specialized  algorithm  for  the  column  update 
can  be  developed  from  an  algorithm  to  efficiently  solve  the  system  Bqyq  =  ce» 
where  <  Vq.Eq  >  is  either  a  rooted  tree  or  a  one-tree. 

If  <  Vq,Eq  >  is  a  rooted  tree,  then  Bq  is  triangularizable  and  yq  can  be 
found  by  forward  substitution.  Consider  the  example  illustrated  in  Figure  18 
with 

1  0  0  0  J136  0 

0-100  _  5 

0  -1  4  0  1/34  0 

2  0  2  1  J  L  J/30  J  L°. 

From  the  first  equation,  5/3$  =  0.  From  the  second  equation,  5/47  =  -5.  From 
the  third  equation,  j/34  =  —5/4  and  from  the  fourth  equation,  j/m  ss  5/2. 
Consider  the  one-tree  illustrated  in  Figure  19  with 
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Figure  about  here 

Figure  18:  Column  update  for  a  rooted  tree  component  of  a  basis  for  a  gener¬ 
alized  network 


The  system  of  equations  for  a  cycle  is  almost  lower  triangular  (only  the  last 
column  has  a  nonzero  entry  above  the  diagonal)  and  takes  the  special  form: 


'  Ol 

fc„  ' 

bi 

a2 

&2 

03 

&3 

fln-1 

_ 

^n  — 1 

a„  _ 

y  =  aei 


which  corresponds  to  the  cycle  illustrated  in  Figure  17.  Under  the  assumption 
that  the  above  system  has  full  row  rank,  the  unique  solution  is  given  by: 


and 


!/i 


ctja\ 

1  -  wiw2-  ■  u„ 


y.+i 


for  i  —  1 , . . . ,  n  —  1 , 

“i+l 


(5) 

(6) 


where  it,  =  -bi/a,  for  i  =  1, . . . ,  n  . 

For  the  example  illustrated  in  Figure  19,  wi  =  l,ti’2  =  1/2, u’3  =  -l,yj2  = 
10/3,3/25  =  -5/3,  and  j/15  =  5/6. 


Figure  about  here 

Figure  19:  Column  update  for  a  one-tree  component  of  a  basis  for  a  generalized 
network 


An  extension  of  the  procedure  cycle  trace  given  in  Section  3.3  can  be  de¬ 
veloped  for  the  case  of  generalized  networks.  The  basic  idea  is  that  the  rooted 
trees  are  triangularizable  and  the  one-trees  are  nearly  triangularizable.  For  the 
components  corresponding  to  one-trees,  (5)  and  (6)  are  used  for  calculations 
involving  a  lower  triangular  system  of  equations.  If  the  nonbasic  column  of  in¬ 
terest,  G  j,  has  two  nonzero  entries,  then  the  procedure  is  applied  twice  and  the 
results  are  added. 


4.4  A  Data  Structure 

The  data  structure  presented  in  Section  3.4  has  been  extended  to  accomodate 
the  one-tree  components  of  a  generalized  network  basis.  With  respect  to  the 


distance  label,  the  nodes  in  the  cycle  of  a  one-tree  and  the  root  node  in  a  rooted 
tree  are  treated  in  the  same  way  with  all  nodes  in  a  cycle  having  distance  label 
of  zero.  The  predecessor  label  of  nodes  in  a  cycle  point  to  the  next  node  in  the 
cycle.  Beginning  with  any  node  in  a  cycle,  say  v,  the  sequence  v,p(v),p(p(v)), . . . 
will  identify  all  nodes  in  the  cycle  and  will  eventually  return  to  node  v.  For  root 
nodes,  v  we  adopt  the  convention  p( v)  =  v.  The  thread  has  also  been  extended 
in  an  obvious  way  with  the  convention  that  traversal  around  a  cycle  using  the 
thread  is  in  the  opposite  direction  to  that  using  the  predecessor.  The  data 
structure  corresponding  to  the  generalized  basis  illustrated  in  Figure  20  is  given 
in  Table  3. 


Figure  about  here 


Figure  20:  Basis  for  a  generalized  network  with  the  predecessor  and  thread 
labels  illustrated 


Table  3.  Data  structure  for  the  basis  illustrated  in  Figure 
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Efficient  procedures  have  been  developed  for  updating  this  data  structure 
after  a  basis  exchange.  As  with  the  linear  network  case,  the  flow  or  value  up¬ 
date  and  the  dual  variable  update  is  usually  integrated  with  the  basis  exchange 
update.  All  of  these  specialized  techniques  can  result  in  software  from  one  to 
two  order  of  magnitude  faster  than  general  linear  programming  software.  All 
calculations  involving  the  flows,  dual  variables,  and  updated  columns  must  be 
performed  using  real  arithmetic  as  opposed  to  the  linear  network  case  which 
only  requires  addition  and  subtraction  of  integers  when  the  data  is  integral. 

Some  additional  computational  efficiencies  can  be  achieved  by  scaling  columns 
of  G  so  that  every  column  has  at  least  one  +1  entry.  The  disadvantage  of  using 
this  scaling  is  that  the  resulting  software  cannot  be  easily  incorporated  within 


a  branch-and-bound  algorithmn  which  can  accomodate  integrality  restrictions 
on  some  or  all  of  the  flow  variables. 


5  Multicommodity  Networks 

Multicommodity  networks  arise  in  practice  when  more  than  one  type  of  com¬ 
modity  must  share  the  capacity  of  the  arcs  in  a  network.  Typical  examples 
include  Air  Force  cargo  routing  models  in  which  cargo  with  different  origin- 
destination  pairs  must  share  the  capacity  of  the  various  aircraft  which  represent 
arcs  in  the  model.  For  this  model,  a  commodity  represents  all  cargo  with  a 
given  base  of  origin.  The  well-known  LOG  AIR  Model  is  a  model  of  this  type. 

Another  well-known  Air  Force  model  is  the  family  of  Patient  Distribution 
System  (PDS)  Models.  The  PDS  generator  was  developed  to  assist  in  evaluating 
the  airlift  capacity  of  moving  patients  from  a  European  theatre  to  U.S.  hospitals. 
One  of  the  input  parameters  for  the  PDS  generator  is  the  number  of  days  in  the 
model.  A  20  day  model  has  over  100,000  columns  in  the  corresponding  linear 
program. 

Let  G  —<  V,  E  >  be  the  underlying  network  through  which  K  commodities 
will  be  flowing  and  let  the  mxn  matrix  A  denote  the  corresponding  node-arc 
incidence  matrix.  Let  the  m-component  vector  rk  denote  the  requirements 
vector  for  commodity  k  and  the  n-component  vector  xk  denote  the  flows  for 
commodity  k  with  corresponding  unit  costs  and  upper  bounds  of  ck  and  u*. 
respectively.  Let  Yk  =  {a-*  :  Axk  =  rk  ,  6  <  xk  <  uk)  and  assume  that  Yk  J  0. 
Suppose  the  shared  (also  called  mutual)  capacity  for  arc  i  is  6j,  an  entry  of  the 
n-component  vector  b.  The  multicommodity  network  flow  problem  is  to  find  K 
n-component  vectors  x 1,...,xli  such  that 


ckxk  =  mint 


It  is  possible  to  generalize  the  multicommodity  network  model  to  allow  for 
both  commodity  dependent  networks  (Akxk  =  rk)  and  for  multipliers  on  the 
mutual  capacity  constraints  (  Dkxk  <  b  ,  where  each  Dk  is  an  nxn  diagonal 
matrix  with  nonnegative  entries).  It  is  aiso  common  that  the  mutual  capacity 
constraints  do  not  involve  all  n  of  the  arcs.  These  generalizations  present  no 
mathematical  difficulties  in  the  algorithms  to  follow  but  greatly  complicate  the 
notation. 

Consider  the  sample  two  commodity  network  illustrated  in  Figure  21.  The 
matrices  and  vectors  corresponding  to  this  model  are  as  follows: 
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r1  =  [  5  0  0  -5  ] 
r*  =  [4  2  -3  —3  ] 
c1  =  [  1  3  5  7  9  ] 

c2~[  2  4  6  8  0  ] 

u1  =  [  4  2  1  3  3  ] 

u2  =  {  2  3  2  2  3  ] 

6= [ 6  3  6  5  3 ] 


Figure  about  here 


Figure  21:  Example  multicommodity  (two  commodity)  network  model 

After  adding  slack  variables  to  the  mutual  capacity  constraints  the  system 
of  equations  corresponding  to  this  model  is 


where  M  is  the  matrix 
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In  general,  the  linear  programming  constraint  matrix  corresponding  to  a 
multicommodity  network  model  takes  the  form: 


Hence,  the  columns  either  take  the  form  for  some  i  or  have  exactly  three 
nonzero  entries  (two  +ls  and  one  —1). 

There  are  three  types  of  specialized  algorithms  which  have  been  used  to 
exploit  the  underlying  structure  of  the  multicommodity  model,  primal  parti¬ 
tioning,  pnce-dtreciive  decommposiiion.  and  resource-directive  decomposition. 
Primal  partitioning  is  a  specialization  of  the  primal  simplex  method  which  ex¬ 
ploits  the  network  structure  of  the  basis.  Price-directive  decomposition  is  a 
specialization  of  Dantzig- Wolfe  decomposition  which  exploits  the  network  struc¬ 
ture  of  the  subproblems.  Resource-directive  decomposition  allocates  the  mutual 
capacity  among  the  A  commodities  and  solves  K  linear  network  problems  per 
allocation.  The  trick  is  to  develop  a  reallocation  scheme  which  will  guarantee 
that  an  optimum  can  be  found. 


5.1  Primal  Partitioning 

By  appending  a  root  arc  to  the  network  for  each  of  the  K  commodities,  we 
convert  the  matrix  (  7)  into  a  matrix  having  full  row  rank.  That  is.  the  fuil 
rank  linear  programming  constraint  matrix  takes  the  form: 


n  +  1  n  +  1 
m  {  A\eP 
m  A\(p 


n  +  1  n 


n  \  /|0  /|0  •••  / 10  // 

where  A  is  a  node-arc  incidence  matrix.  It  is  well  known  that  every  basis  for 
(  8)  may  be  partitioned  in  the  following  form: 
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where  B1, _ Bh  correspond  to  rooted  trees.  This  basis  takes  the  general  form 


where 


m  A  n  —  p  p 

mK  /  L\  R\ 

q  \  L-2  R2 

P  \  L3  R3  1 


R7-[TX 

l3=[s1 


*3=[  U'  ...  V«  ) 

By  partitioning  the  dual  variables  to  be  compatible  with  (  9),  the  dual  calcula¬ 
tion  requires  the  solution  to  the  following  system: 


rrj  |  7t2  |  ff3 


Lx 

Rx 
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l2 

r2 

Lz 

Rz 

~r 

=  (  ci  c2  0  ] 


Let  Q  —  [ R2-L2Ll  1  /2i].  Then  tr3  =  6,  *2  =  (c2-ctL~'  Ri  )Q~l,  and  Jr2  =  (cj  - 
5 t2L2)LJ"1.  Note  that  L\  corresponds  to  K  rooted  trees  so  that  all  calculations 
involving  LJ"1  can  be  executed  directly  on  the  rooted  trees  and  do  not  require 
the  use  of  explicit  matrices.  The  only  matrix  operations  involve  the  qxq  matrix 
Q  which  is  called  the  working  basts. 

Partitioning  the  updated  column  in  a  similar  way,  the  column  update  calcu¬ 
lation  requires  the  solution  of  the  following  system: 


Li  Ri 
L2  r2 
l3  r3  I 
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The  solution  is  given  by  y  =  <2-1(r  ~  L^L^u),  x  —  Lj"!u  -  If1  R\y,  and 
z  =  tt>  -  L3L^lu  -  (R3  -  LsL^R^y-  As  before,  the  only  matrix  operations 
ir  .olve  the  working  basis  Q. 

Each  basis  exchange  can  alter  the  dimension  of  the  working  basis  by  at  most 
one.  Efficient  techniques  have  been  developed  for  maintaining  Q~1  as  the  basis 
is  updated.  If  q  is  small,  then  most  of  the  work  associated  with  a  simplex 
pivot  can  be  performed  efficiently  using  the  K  rooted  trees.  If  q  is  large,  then 
the  pivots  can  be  more  expensive  than  a  regular  simplex  pivot  which  does  not 
exploit  the  special  network  structure. 

5.2  Price-Directive  Decomposition 

For  each  k£tK ,  suppose  that  Yk  =  :  Axk  =  rk  ,  6  <  zk  <  u*}  ^  0 

and  is  bounded,  and  let  tlf, . . .  ,u,kt  denote  the  extreme  points  of  Yk .  Let  the 
nxn*  matrix  IF*  be  given  by  [  wk  j  wk  \  ■  ■  ■  j  u>,*t  j .  Then  any  xk6Yk  can  be 
represented  as  a  convex  combination  of  those  extreme  points  (the  columns  of 
\Vk).  That  is,  xk  =  IF*  A*  where  1A*  =  1  and  A*  >  0.  Let 

W  =  {(  A1  |  ••  |  A*  j  s  )  :  ]jTtW'iA*-M  =  f>,  s  >  6 ,  1  A*  =  1 ,  A*  >  6} 

Then  an  equivalent  statement  of  the  multicommodity  network  problem  is  to 
find  A'  +  1  vectors  (  A1  |  • • •  |  A*  |  s  )  such  that 

Y)ickWlXi  =  min{J^kckWkXk  :  (  A1  |  •  •  |  A*  |  s  )eIF} 

The  disadvantage  of  this  statement  of  the  problem  is  that  it  is  quite  difficult  to 
find  all  the  extreme  points  for  Yk  (which  form  the  columns  of  IF*).  However,  the 
price-directive  algorithm  provides  a  mechanism  for  applying  the  primal  simplex 
algorithm  to  this  formulation  by  providing  columns  of  VF*  only  when  they  are 
needed.  In  fact,  the  complete  matrix  IF*  is  never  actually  generated. 

Suppose  we  have  a  nonnegative  basic  feasible  solution  to  the  system  of  equa¬ 
tions: 


^  IF*A*  +  s  =  6  (x) 

k 

iA1  =  1  (7.) 

ixK  =  1  (7k) 

with  the  corresponding  duals  (  tr  |  71  |  ■  ■  ■  |  7*  )•  From  Section  2,  we  see  that 
the  tth  nonbasic  slack  variable  prices  favorably  if  >  0.  Also,  A*  prices  favor¬ 
ably  if  wkx  +  7k -ckuj  >  0.  Note  that  the  extreme  point  from  Yk  which  prices 
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most  favorably  can  be  found  by  solving  the  linear  network  problem 

min{(ci  -  x)u*  :  wk^Yk} 

using  an  extreme  point  algorithm.  If  the  extreme  point  produced  from  solving 
the  linear  network  problem  prices  unfavorably,  then  none  of  the  columns  of 
Wk  will  price  favorably.  Hence,  columns  of  Wk  are  generated  only  as  they  are 
needed. 

The  linear  program 

minjJ^W^A*  :  (  A1  |  •  •  •  |  \k  f  s  )  €W) 

is  known  as  the  master  problem  and  the  linear  network  problems  minfc**1  : 
z*€V*}  for  k&K  are  called  the  slave  subproblems.  Columns  of  the  master  prob¬ 
lem  (extreme  points  of  the  V'*)  are  found  by  solving  linear  network  problems. 
A  mathematical  description  of  the  price-directive  algorithm  follows: 


procedure  PRICE-DIRECTIVE  DECOMPOSITION 


inputs:  A 

mxn  node-arc  incidence  matrix 

K 

-  number  of  commodities 

cK... 

,cK 

n-component  vectors  of  unit  costs 

uK... 

n-component  vectors  of  upper  bounds 

r1,... 

,r“ 

-  m-component  vectors  of  requirements 

6 

n-component  vector  of  mutual  capacities 

outputs:  y1,... 

,yK 

n-component  vectors  of  optimal  flows 

assumption:  Yk  = 

. 

Axk  =  r*,0  <  xk  <  u*}  ^  $  and  is  bounded 

for  all  k  £  tK. 

begin 

[ Initialize  Master  Problem] 

obtain  a  basic  feasible  solution  [y1 1  - .  •  | yK]  for  W  with  corresponding  dual  variables 
Ml]  (if  a  basic  feasible  solution  is  not  readily  available,  then  artificial  variables  and 
a  two-phase  procedure  may  be  used); 
optimal  <=  ‘no’; 
while  optimal  =  ‘no’  do 
[Price  Slacks] 

t  <=  1; 

while  i  <  n  do 
if  *  <  <  0 
then 

perform  one  simplex  iteration  in  the  master  problem  using  [e,  |0]  as 
the  nonbasic  column  which  prices  favorably; 
update  [yl|.  -li/*]  and  [fffr]; 

*  ■<=  1; 
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else 

t  <=  i  +  1; 
endif 
end  while 

[Price  Extreme  Points] 
k  <=  1,  favorable  ‘no’; 
while  k  <  K  and  favorable  =  ‘no’do 
obtain  an  optimal  extreme  point  wk  for  min{(c*  -  *■)**  :  Axk  = 
rk ,  0  <  xk  <  u*}; 
if  (c*  -  x)xk  <  7* 
then 

favorable  <=  ‘yes’; 
else 

k  <S=  k  4-  1; 
endif 
endwhile 
if  favorable  =  ‘yes' 
then 

perform  one  primal  simplex  iteration  in  the  master  problem  using 
[tc*|e*]  as  the  nonbasic  column  which  prices  favorably; 
update  [y1 1  -  - .  |yx]  and  [x|-y]; 
else 

optimal  =  ‘yes’; 
endif 
endwhile 
end 

A  pictorial  description  of  the  price-directive  decomposition  algorithm  is  given 
in  Figure  22.  The  master  program  produces  dual  variables  used  to  construct  the 
objective  function  coefficients  for  each  of  the  slave  subproblems.  Each  slave  sub- 
problem,  in  turn,  produces  an  extreme  point  which  may  or  may  not  provide  an 
objective-improving  column  relative  to  the  current  basis  for  the  master  problem. 
A  favorable  column  will  be  used  for  a  simplex  iteration  in  the  master  problem, 
producing  new  dual  variables  for  the  slave  subproblems.  Note  that  successive 
calls  to  a  routine  which  solves  a  given  subproblem  involve  only  a  change  in  the 
objective  function  and  do  not  involve  any  changes  to  the  constraints.  Hence 
each  such  call  can  make  use  of  the  optimal  basis  tree  produced  by  the  previous 
call.  We  also  recommend  that  a  wrap-around  list  be  maintained  to  organize  the 
calls  to  the  subproblems,  a  minor  variation  on  the  above  procedure.  Tha*  «s,  if 
commodity  k  produced  a  favorably  priced  column  at  iteration  i,  then  begin  the 
new  iteration  by  pricing  the  columns  associated  with  commodity  (k+  l)modA’. 
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Figure  about  here 


Figure  22:  The  price-directive  decomposition  algorithm 


5.3  Resource-Directive  Decomposition 

It  is  the  mutual  capacity  constraints  xk  <  b  that  make  the  multicommodity 
problem  difficult.  If  these  constraints  could  be  eliminated  or  ignored,  then  the 
problem  would  decompose  into  K  linear  network  problems.  Let  [  y1  \  ■  ■  ■  \  y*  ] 
denote  an  allocation  of  b  among  the  K  commodities  such  that  0  <  y*  <  u*  for 
each  k  and  y*  <  b.  Let 

zl( yk)  =  min{clxk  :  Axk  =  rk  ,  0  <  xk  <  y*} 

An  equivalent  statement  of  the  multicommodity  network  problem  is  to  find  K 
vectors  [  y1  |  • ••  |  y*  ]  such  that 

53-*(y*)  =  min{^  ckxk  :  Axk  =  rk  ,  6  <  xk  <  yk) 
k  k 

The  objective  function  y([  y1  |  •  •  •  |  y*  ] )  =  zk(yk)  is  piece-wise  linear  and 
convex.  Hence,  the  classic  subgradient  algorithm  can  be  applied  to  this  problem. 

Consider  the  nonlinear  program  min{y(y)  :  yGV),  where  g(  )  is  piece-wise 
linear  and  convex  ai.d  Y  £  9  is  formed  by  the  intersection  of  a  finite  number 
of  closed  half-spaces.  To  apply  the  subgradient  algorithm,  one  must  be  able  to 
solve  the  problem 

ntinfQTfy,  -p<)2]l/2  ;  y€Y] 

i 

for  any  point  p.  The  solution  of  this  problem  is  called  the  projection  of  p  onto 
Y.  This  is  traditionally  expressed  as  y  < —  P\p]- 

Mathematically  the  subgradient  algorithm  may  be  stated  as  follows: 


procedure  SL’BGRADIENT 


inputs: 


output: 

assumptions: 


g(  )  -  the  objective  function 

y  -  the  feasible  region 

a,,s 2,  --  -  the  sequence  of  step  size  values 

y  -  an  optimum  for  min{y(y)  :  y  €  Y) 

y(  )  is  piece- wise  linear  and  convex 

y  5*  <J>  is  formed  by  the  intersection  of  a  finite  number  of 
closed  half-spaces 


begin 

[Initialization] 

obtain  a  point  y  6  y,  optimal  <=  ‘no’,  i  <=  1; 


while  optimal  =  ‘no-  do 
obtain  a  subgradient  ti  of  g(  )  at  y; 
if  ti  ^  6 
then 

optimal  <=  ‘yes’; 
else 

y  <=  P[y  -  » <= » + 1; 

endif 

endwhile 

end 

Several  convergence  results  for  this  algorithm  have  appeared  in  the  literature. 
These  results  all  provide  restrictions  and  guidance  on  the  selection  of  the  se¬ 
quence  of  step  sizes. 

The  computationally  expensive  part  of  the  above  algorithm  when  applied  to 
the  multicommodity  network  flow  problem  is  the  calculation  of  the  subgradient. 
That  is,  we  need  a  subgradient  of  the  function  g(  [  y 1  |  •  •  •  |  y*  ] )  at  the  point 
[  yl  |  •  •  •  |  i/K  ]  ,  where  g([  yl  |  •  •  j  yK  ])  =  **($*)•  By  duality  theory, 

zk(yk)  =  min{ctii  :  Azk  =  rk  ,  6  <  xk  <  yk) 

=  ma x{rkyk  -  ykvk  :  ykA  -  i/k  <  ck  ,  i/k  >  6} 

Let  for  all  i€*A"  denote  optimd  dual  variables  for  ff([  y1  |  •  ••  |  yK  ])- 

It  can  easily  be  shown  that  [  -P1  |  •  •  •  |  —PK  ]  is  a  subgradient  for  g(  )  at  the 
point  ( [  y1  |  •  •  •  |  yK  ] )  .  Efficient  projection  algorithms  for  the  special  case  of 
the  multicommodity  network  problem  are  also  available  in  the  literature. 

The  major  attraction  of  the  subgradient  algorithm  is  that  the  multicommod¬ 
ity  network  problem  can  be  solved  by  solving  a  sequence  of  single  commodity 
network  problems.  The  major  disadvantage  is  that  it  is  quite  difficult  to  select  a 
set  of  step  sizes  so  that  the  software  implementation  is  robust  over  a  wide  range 
of  problems.  Our  experience  has  indicated  that  much  skill  is  required  in  the 
selection  of  the  step  sizes.  Many  skilled  mathematical  programming  software 
developers  have  discovered  that  the  step  size  selection  is  more  of  an  art  than  a 
science. 


6  Reference  Notes 

The  specialized  algorithms  presented  in  this  chapter  all  rely  on  the  graphical 
interpretation  of  the  simplex  steps  when  applied  to  a  linear  program  possessing 
an  underlying  network  structure,  in  whole  or  in  part.  Graphical  characteriza¬ 
tions  for  bases  for  both  the  linear  network  problem  and  the  generalized  network 
problem  can  be  found  in  Dantzig  [1963].  Additional  elaboration  on  the  inter¬ 
pretation  of  the  algebraic  operations  on  a  graphical  structure  were  provided  by 


Johnson  [1965].  The  first  software  implementations  which  empirically  demon¬ 
strated  the  merit  of  those  specialized  algorithms  were  developed  by  Srinivasan 
and  Thompson  [1973]  ar.  1  by  Glover  and  Klingman  and  their  colleagues  at  the 
University  of  Texas  (see  Glover.  Karney,  and  Klingman  [1974],  Glover,  Karney, 
Klingman,  and  Napier  [1974],  and  Glover,  Klingman,  and  Stutz  [1974]  ).  Since 
these  seminal  papers,  hundreds  of  papers  with  various  extensions  and  special¬ 
izations  have  appeared  in  the  literature.  Bradley,  Brown,  and  Graves  [1977]  and 
Barr,  Glover,  and  Klingman  [1979]  are  two  of  our  favorites. 

The  first  specialized  primal  simplex  code  which  exploited  the  graphical  na¬ 
ture  of  the  basis  of  a  generalized  network  was  developed  by  Glover,  Klingman, 
and  Stutz  [1973].  Since  then,  many  excellent  software  implementations  have 
been  developed.  A  table  of  codes  and  references  may  be  found  in  Clark,  Ken- 
nington,  Meyer,  and  Ramamurti  [1992].  One  of  our  favorite  generalized  network 
codes  is  GENNET  which  was  developed  by  Brown  and  McBride  [1984]. 

The  primal  partitioning  method  for  multicommodity  problems  is  due  pri¬ 
marily  to  Hartman  and  Lasdon  [1970,1972].  The  price-directive  decomposition 
procedure  was  first  developed  by  Ford  and  Fulkerson  [1958].  It  is  purported 
that  the  well-known  Dantzig- Wolfe  decomposition  algorithm  (see  Dantzig  and 
Wolfe  [1960]  )  was  inspired  by  the  Ford  and  Fulkerson  paper  on  multicommodity 
problems.  The  subgradient  algorithm  was  first  applied  to  the  multicommodity 
network  problem  by  Held,  Wolfe,  and  Crowder  [1974]. 

In  this  chapter,  the  algorithms  presented  were  generally  based  on  a  special¬ 
ization  of  the  simplex  method.  There  are  three  other  basic  approaches  which 
are  now  competing  quite  successfully  in  empirical  analyses.  The  oldest  is  the 
relaxation  method  which  was  developed  by  Bertsekas  and  his  colleagues  at  MIT 
(see  Bertsekas  [1991]  ).  The  second  is  the  algorithm  of  Goldberg  which  has  been 
implemented  in  software  by  Anderson  and  Setubal  [1992].  The  third  algorithm 
which  could  have  an  impact  on  the  field  is  the  network  interior  point  method 
of  Resende  and  Veiga  [1992].  All  of  these  approaches  have  their  advocates  and 
we  expect  improvements  for  each  of  these  relatively  new  algorithms  to  appear 
in  the  literature  in  the  near  future. 

Many  books  exist  which  contain  excellent  presentations  of  the  various  algo¬ 
rithms  for  a  variety  of  network  problems.  Our  favorites  include  Papadimitriou 
and  Steiglitz  [1982],  Chvatal  [1983],  and  Ahuja,  Manganti,  and  Orlin  [1989]. 
This  chapter  would  not  be  complete  without  mentioning  our  own  book  (Ken- 
nington  and  Helgason  [1980])  which  contains  a  wealth  of  technical  details  not 
found  in  other  publications. 
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(a)  Generated  subgraph  G({1,2,3}) 


Figure  2:  Example  network  subgraphs 


P(a)  =  {4, C43.3, €*,4,624, 2} 

(a)  Nonpath  walk  from  node  4  to  node  2 


PW  =  {2,€23,3,e34,4,€24i2} 


(c)  A  cycle  on  nodes  {2,3,4} 


Figure  3:  Walks  in  the  example  network 
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(b)  Basis  tree  2 

Figure  7:  Basis  trees  for  sample  network 
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(a)  Costs  for  sample  basis  tree 


(b)  Dual  variables  for  sample  basis  tree 


Figure  8:  Dual  variable  calculation  example 
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(a)  Simple  path  linking  5  to  6 
:  yA  =  C 


Figure  10.  Updated  column  example 
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Figure  16:  One-tree  component  of  a  basis 
for  a  generalized  network 
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Figure  22:  The  price-directive  decomposition  algorithm 


APPENDIX  B 


The  One-to-One  Shortest-Path  Problem:  An  Empirical 
Analysis  With  the  Two- Tree  Dijkstra  Algorithm 

Richard  V.  Helgason  t 
Jeffery  L.  Kennington  f 
B.  Douglas  Stewart  t 

Four  new  shortest-path  algorithms,  two  sequential  and  two  parallel,  for  the  source  to  sink  shortest- 
path  problem  are  presented  and  empirically  compared  with  five  algorithms  previously  discussed 
in  the  literature.  The  new  algorithm,  S22,  combines  the  highly  effective  data  structure  of  the 
S2  algorithm  of  Dial,  Glover,  Karney,  and  Klingman,  with  the  idea  of  simultaneously  building 
shortest-path  trees  from  both  source  and  sink  nodes,  and  was  found  to  be  the  fastest  sequential 
shortest-path  algorithm.  The  new  parallel  algorithm,  PS22,  is  based  on  S22  and  is  the  best  of  the 
parallel  algorithms.  We  also  present  results  for  three  new  S22-type  shortest-path  heuristics.  These 
heuristics  find  very  good  (often  optimal)  paths  much  faster  than  the  best  shortest-path  algorithm. 

Revised  (December  1992) 


Since  the  late  fifties  when  its  first  solution  methods  were  developed,  the  shortest-path  problem  has 
become  one  of  the  fundamental  problems  in  the  areas  of  combinatorial  optimization,  computer  science,  and 
operations  research.  Algorithms  and  applications  are  commonly  found  in  the  important  books  in  these  areas 
(see  for  example  Berge  and  Ghouila-Houri  (1962),  Bertsekas  and  Gallager  (1987),  Even  (1979),  Hu  (1982), 
Jensen  and  Barnes  (1980),  Lawler  (1987),  Papadimitriou  and  Steiglitz  (1987),  Quinn  (1984),  Rockafellar 
(1984),  and  Tarjan  (1983)).  The  study  of  this  problem  has  been  motivated  by  both  its  elegant  mathematical 
structure  and  its  many  practical  applications.  Our  recent  interest  in  this  problem  was  occasioned  by  the 
need  to  solve  shortest-path  subproblems  in  several  mathematical  optimization  procedures  we  are  developing 
in  an  MIMD  parallel  computing  environment. 

Excellent  surveys  of  the  many  shortest-path  problem  variations  may  be  found  in  Deo  and  Pang  (1984) 
and  Gallo  and  Pallottino  (1986).  A  survey  of  techniques  and  computational  comparisons  may  be  found  in 
Gallo  and  Pallottino  (1988),  in  Dial,  Glover,  Karney,  and  Klingman  (1977)  and  (1979),  in  Klingman,  Mote, 
and  Whitman  (1978),  in  Glover,  Glover,  and  Klingman  (1984),  in  Desrochers  (1987),  and  in  Divoky  (1987). 
The  methods  are  grouped  into  two  general  classes:  label-setting  algorithms  and  label-correcting  algorithms. 
Dijkstra  (1959)  is  credited  with  the  first  label-setting  algorithm  and  any  algorithm  that  uses  this  approach 
has  been  considered  a  particular  implementation  of  Dijkstra ’s  original  algorithm  (see  Gallo  and  Pallottino 
(1986)). 

Typically,  the  shortest-path  problem  is  one  that  requires  the  shortest-path  from  a  single  source  node,  say 

s,  to  all  other  nodes  in  a  network.  The  solution  can  be  represented  as  a  shortest-path  tree  rooted  at  s.  In  this 
paper  we  are  concerned  with  the  problem  of  finding  the  shortest-path  between  a  source  node  and  a  sink  node, 

t.  Dantzig  (1960)  suggested  that  a  pair  of  trees  be  built  with  one  rooted  at  s  and  the  other  rooted  at  t.  No 
stopping  criteria  were  given.  This  strategy  also  appears  in  the  book  by  Berge  and  Ghouila-Houri  (1962)  with 
an  incorrect  stopping  criterion.  Nicholson  (1966)  was  the  first  to  present  a  correct  analysis  of  the  Dijkstra 
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two-tree  algorithm.  Hart,  Nilsson,  and  Raphael  (1968)  presented  a  one-tree  algorithm  utilizing  heuristic 
cost  functions.  They  included  the  case  in  which  lower  bounds  on  distances  between  nodes  are  available  and 
proved  that  optimal  paths  are  obtained.  Pohl  (1971)  extended  these  results  for  the  bi-directional  algorithm, 
antedating  Mohr  and  Pasche  (1988),  who  presented  similar  results.  Additional  discussion  may  be  found  in 
the  survey  by  Dreyfus  (1969),  where  he  conjectured  that  building  trees  from  s  and  t  would  be  ineffective  for 
this  problem. 

Few  empirical  studies  have  been  reported  in  the  literature  for  the  two-tree  algorithms.  Pohl  (1969) 
presented  results  that  have  received  little  recognition,  and  no  further  studies  appear  until  recently.  Mohr 
and  Pasche  (1988)  solved  for  shortest-paths  in  three  grid  networks  and  a  network  representing  a  road  map 
of  Switzerland,  using  the  two-tree  Dijkstra  algorithm  as  well  as  their  version  of  the  two-tree  algorithm  for 
networks  with  lower  bounds  available.  They  also  simulated  results  for  parallel  algorithms  if  two  processors 
were  available. 

The  motivation  for  this  study  was  to  implement  two-tree  algorithms  using  more  efficient  data  structures 
than  the  classical  Dijkstra  algorithm,  and  to  actually  solve  shortest-path  problems  using  two  processors. 
Four  new  algorithms  were  developed:  sequential  and  parallel  two-tree  algorithms  based  on  Dial  (1969)  (the 
Si  code  in  Dial  et  al.  (1979)),  and  sequential  and  parallel  two-tree  algorithms  based  on  the  S2  code  in  Dial 
et  al.  (1979).  These  codes  are  compared  with  the  classical  Dijkstra,  SI,  S2,  two-tree  Dijkstra,  and  parallel 
two-tree  Dijkstra  codes. 

While  performing  this  study,  it  was  noted  that  the  two-tree  S2  algorithm  finds  “good”  paths  quickly. 
Three  heuristic  path  algorithms  were  developed  by  simply  stopping  early  while  executing  the  two-tree  S2 
algorithm.  Often  times  these  heuristics  find  a  shortest-pat-  ,  although  they  do  not  prove  it  is  so,  and  these 
paths  are  found  much  faster  than  the  fastest  optimal  algorithm.  We  present  computation  times  as  well  as 
measures  of  closeness  to  optimality  for  these  heuristics. 

1.  THE  ALGORITHMS 

This  section  presents  the  definitions,  notation,  and  nine  algorithms  representing  the  codes  used  in  our 
computational  study.  The  notation  and  presentation  are  based  on  that  found  in  Gallo  and  Pallattino  (1986). 
The  input  for  each  algorithm  is  a  directed  graph  G  =  [/V,A]  with  a  node  set  N  and  an  arc  set  A.  Associated 
with  each  arc  (i,j)  £  A  is  a  length  Uj.  A  shortest-path  is  desired  between  two  nodes  s  and  t  and  it  is 
assumed  that  such  a  path  exists.  We  also  assume  for  all  algorithms  that  >  0  for  all  (*,  j)  £  A.  For 
efficient  implementation  it  is  assumed  that  G  is  given  in  the  form  of  arc-lists.  Specifically,  the  forward  star 
for  node  u  £  N,  FS(u),  is  the  set  of  all  arcs  (u,  j)  £  A.  Six  algorithms  require  a  backward  star,  BS(v ), 
defined  as  the  set  of  all  arcs  (»,  v)  £  A  for  node  v  £  N. 

The  basic  working  entities  of  each  algorithm  include  a  set  of  labels,  du,  for  node  distances  from  the  root 
of  a  shortest-path  tree  T,  a  set  of  predecessors,  p„,  for  nodes  in  the  tree,  and  the  set  of  candidate  nodes, 
Q.  In  the  algorithms  based  on  the  SI  and  S2  algorithms,  the  set  Q  will  be  divided  into  subsets,  or  buckets, 
that  will  contain  candidate  nodes  that  are  the  same  distance  from  the  root  node.  This  requires  that  each 


arc  length  be  mitger  to  correspond  to  an  index.  For  ease  of  presentation,  we  will  let  T  be  the  set  of  indices, 
{0, ... ,lmax }  for  the  buckets  of  Q ,  where  Imax  =  max{/u„  •  {«,  v)  £  A).  The  notation  will  be  modified  lor 
the  bi-directional  algorithms  by  using  the  superscripts  s  or  t  to  indicate  the  root  of  the  tree,  The  algorithms 
terminate  with  the  length  of  the  shortest-path  from  s  to  t  and  a  small  set  of  nodes,  J,  used  to  identify  a 
shortest-path.  A  shortest-path  from  s  to  t  is  implicit  in  the  predecessor  labels. 

1.1  The  classical  Dijkstra  algorithm 

Dijkstra’s  classical  algorithm  (1959)  begins  at  node  s  and  builds  a  shortest-path  tree  T  in  which  the 
shortest-path  from  s  to  any  node  in  the  tree  is  known.  When  node  <  is  placed  in  the  tree  we  have  a  minimum 
length  directed  path  from  s  to  t  and  may  terminate.  As  mentioned  before,  the  algorithms  discussed  here 
differ  in  the  way  the  candidate  nodes  are  placed  in  and  retrieved  from  the  set  Q.  The  implementation  here 
for  the  classical  Dijkstra  algorithm  (Dl)  searches  the  list  of  candidate  nodes,  Q,  for  minimum  label  nodes 
and  places  all  such  nodes  in  the  set  R.  Then  all  the  nodes  in  R  can  be  scanned  one  after  the  other.  When 
there  are  many  nodes  tied  with  the  minimum  label,  searches  are  avoided  with  only  minimal  effort.  The 
algorithm  may  be  stated  as  follows: 

Procedure  Dl(s,f) 
begin 

initialize: 

pi  *-  0 ,di  «—  oo  for  all  i  £  IV;  Q  <—  0; 
d.  —  0,p,  «-  s,R  «-  {s},  T  <—  0; 
while  t  £  R 

for  each  u  £  R  do 

for  each  (u,  v )  £  FS(u)  such  that  du  +  lvv  <  dv  do 
dv  <  du  +  luv ; 

Pv  —  «; 

if  v  £  Q  then  Q  <—  Q  U  {t?}; 

end 

r-ru{u}; 

end 

comment:  search  Q  for  minimum  label  nodes  and  place  in  R 
a  *—  min{di  :  i  6  Q},R  <—  {»  :  di=a},Q  Q\R ; 

endwhile 
end 

1.2  The  SI  algorithm 

This  one-tree  algorithm  is  implemented  as  proposed  by  Dial  (1969).  As  in  algorithm  Dl,  it  selects  a 
minimum  label  node  to  scan  each  iteration,  but  the  nodes  in  Q  are  stored  in  buckets  according  to  their 
distance  labels.  More  specifically,  Imax  +  1  buckets  are  required,  where  Imax  =  max{/ov  :  (u,u)  £  A}  and 
a  node  u  is  stored  in  bucket  z  if  du  =  z(mod  Imax  4- 1).  Only  nodes  with  equal  distance  labels  will  be  in 
bucket  z  and  only  Imax  +  1  buckets  are  required,  because  for  a  node  u  with  minimum  label,  we  have  that 
for  each  v  €  Q,du  <  dv  <  du  +  (max.  The  buckets  are  implemented  efficiently  as  two-way  linked-lists.  From 
minimum  label  node  u,  the  next  non-empty  bucket  contains  the  next  minimum  label  node(s)  and  the  effort 
to  search  an  entire  list  as  in  algorithm  Dl  is  greatly  reduced.  The  trade-off  is  increased  effort  managing 


the  two-way  linked-list  each  time  a  node  has  an  improved  distance  label.  The  algorithm  may  be  stated  as 
follows: 

Procedure  Sl(a,i) 
begin 

initialize: 

Pi  <—  0,d<  ♦—  oo  for  all  i  e  N‘,  Qj  ♦—  0  for  z  =  l,...,imaz; 

Qo  —  {*}.  d,  «-  0 ,p,  ♦-  s,  T  ♦—  0; 
while  t  g  T 

let  z  be  the  next  index  such  that  Q,  ^  0  ; 
for  each  u  €  Q,  do 

Q,  -  9,\W; 

for  each  (u,  v)  6  FS(u)  such  that  du  +lvv  <  do 
a  «—  dy  (mod  Imax  +  1); 
dv  *—  du  +  luv ; 

6  <—  <fv(mod  Imax  +  1); 
pu  ♦-  u; 

Qa  QoM*'}; 

Qb  Qi  u  {v}; 

end 

T «—  TU  {«}; 

end 

endwhile 

end 


1.3  The  S2  algorithm 

This  one- tree  algorithm  is  based  on  an  idea  due  to  Dantzig  (1960)  and  the  implementation  here  is  due 
to  Dial  et  al.  (1977).  It  requires  that  each  FS(u)  for  all  u  €  N  be  sorted  in  shortest  first  order.  Given  this, 
the  observation  can  be  made  that  the  entire  forward  star  of  ?  Tied--  u  need  not  be  scanned  all  at  once  in  that 
the  node,  say  v,  that  is  first  updated  will  have  a  distance  label  less  than  or  equal  to  any  subsequent  nodes 
updated  from  node  u.  The  node  v  is  placed  on  a  one-way  linked-list,  paired  with  u,  at  a  level  du  +  /u„.  In 
stating  the  algorithm  below,  we  use  ifc(u)  as  a  counter  to  point  to  the  next  arc  in  FS(u)  to  be  scanned.  The 
node  u/*(u)  is  the  corresponding  node  of  this  arc.  The  length  of  each  forward  arc-list  is  given  by  A(u).  It 
should  be  noted  that  the  actual  implementation  does  not  use  a  counter.  A  simple  minus  sign  in  the  forward 
star  array  can  indicate  the  next  arc  to  be  scanned. 

Nodes  may  be  duplicated  on  the  linked-list,  therefore  no  deleting  is  required.  Even  so,  the  number  of 
nodes  on  the  linked-list  will  not  exceed  the  total  number  of  nodes  since  only  one  per  scanned  node  is  allowed. 
Q  is  searched  as  in  algorithm  Si  for  the  next  non-empty  bucket  and  minimum  label  node.  If  the  node’s 
paired  predecessor  is  not  its  current  predecessor,  the  node  is  already  in  the  shortest-path  tree  and  may  be 
discarded.  The  algorithm  may  be  stated  as  follows  (see  Gallo  and  Pallottino  (1986)): 


Procedure  S2(s,t) 

begin 

initialize: 

Pi  <—  0 ,di  «—  oo,k(i)  *—  1,  A(i)  =  |F5(»)|  for  all  i  €  N; 

Q,  *—  0  for  z  =  1, Imai;  Qo  *-  {s},  d,  *- Q,p,  —  a,  T  *—  0; 
while  <  £  T 

let  z  be  the  next  index  such  that  Q,  £  0; 
for  each  u  6  Q,  do 

comment:  determine  next  node  in  FS(u) 

V  <-  Wk(u)\ 

Q *  —  QAM; 

comment:  determine  new  candidate  arc 

INSERTM; 

r-TU{tij; 

if  t>  i  Q  then  INSERT^); 
end 

endwhile 

end 

Procedure  INSERT(x) 
begin 

Jfc(z)  «—  k(x)  +  1; 

y  *-  «<*(r); 

while  k(x)  <  h(x)  and  dx  +  lxy  >  dy  do 
k(x)  <-  k(x)  +  1; 

y «- 

endwhile 
if  k(x)  <  h(x)  then 
Py  —  x; 

dy  dx  +  fry'. 

a  <—  dy(mod  Imax  -f  1); 

Qa  —  C?aU{z}; 
endif 
end 


1.4  The  two- tree  Dijkstra  algorithm 

The  two-tree  Dijkstra  algorithm  builds  a  pair  of  shortest-path  trees,  one  from  s  and  one  from  t.  The 
tree  rooted  at  t  is  analogous  to  the  one  rooted  at  s,  but  scans  backward  stars,  and  its  predecessors  are  the 
heads  of  arcs  rather  than  tails.  The  two  trees  are  grown  in  alternate  steps  and  termination  is  triggered 
when  a  node  appears  in  both  trees.  It  should  be  noted  however  that  the  node  appearing  in  both  trees  is  not 
necessarily  on  a  shortest-path  (see  Nicholson  (1966)  for  a  proof  of  termination  criteria).  Testing  on  random 
graphs  showed  that  about  72%  of  the  time  this  node  will  be  on  the  shortest-path.  A  search  is  performed 
over  nodes  in  both  trees  to  find  a  node,  say  r,  such  that  d*  +  d*r  is  a  minimum.  A  shortest-path  from  s  to  t 
may  be  found  by  following  predecessors  from  r  to  s  and  from  r  to  t  in  each  tree,  respectively.  We  define  J 
as  the  set  of  nodes  which  can  be  used  to  identify  a  shortest-path.  The  algorithm  requires  twice  the  storage 
of  algorithm  D1  and  may  be  stated  as  follows: 
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Procedure  D2(s,<) 
begin 

initialize: 

Pi  «—  0,d{  —  oo  for  all  i  6  A'';  Q*  <—  0; 
dJ-0 

p\  *—  0,  df  *—  oo  for  ail  i  €  N;  Q*  —  0; 

4  —  0,pj  t,r  -  0,  R*  4-  {t}; 

while  t* nr1  =  0  do 
for  each  u  €  R*  do 

for  each  (u,  v)  6  FS(u)  such  that  d*u  +  <  dj  do 

dj  dl  +  luv; 
pi  — 

if  t>  £  Q*  then  Q*  <—  (J*  U  {v}; 
end 

T'  *—T*  U  {u}; 
end 

comment:  search  Q‘  for  minimum  label  nodes  and  place  in  R* 
a  -  min{<  :  i  6  Q F*  4-  {,' :  d\=a},Q*  - 

for  each  v  E  R*  do 

for  each  (u,  t>)  €  BS(u)  such  that  dj,  +  <  d*,  do 

dj,  dj,  +  Iuv ; 

pi  v; 

if  v  £  Q*  then  Q*  *-  Q*  U  {u}; 
end 

T*  4—  T*  U  {«}; 
end 

comment:  search  Q*  for  minimum  label  nodes  and  place  in  R ‘ 
«  -  minK  : »  €  Q‘},  R f  4-  {t  :  Q*  <-  Q‘\ie‘; 

endwhile 

comment:  stopping  criterion  met 
f3  4—  min{dj  +  d}  :  i  ET‘  UT*}; 

7  4—  {i  €  T*  UT*  :  dj  +  d|=/?}; 
end 


1.5  The  two-tree  Si  algorithm 

The  two-tree  SI  algorithm  builds  a  pair  of  shortest-path  trees,  one  from  s  and  one  from  t,  using  SI  data 
structures  for  each  tree.  As  in  algorithm  D2,  termination  is  triggered  when  a  node  is  first  found  to  be  in 
both  trees.  The  node  r  such  that  d‘r  +  dj.  is  minimum  gives  the  shortest- distance  from  s  to  t  and  lies  on  a 
shortest-path.  The  algorithm  requires  twice  the  storage  of  SI  and  may  be  stated  as  follows: 

Procedure  S12(s,t) 
begin 

initialize: 

Pi  4—  0,  d*  4—  00  for  all  i  6  Q*  *—  0  for  z  =  1, Imax; 

Qo  {*}>  dj  4—  0,pj  4-  s,  T'  4—  0; 

p‘  4—  0,  dj  4—  co  for  all  i  E  N]  4—  0  for  z  =  1, ...,  Imax; 

Qo  -  {<},  d*.  -  O.pJ  -  t,  r  *-  0; 
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while  T*  HT*  =  0  do 

let  z  be  the  next  index  such  that  Q\  £  0  ; 
for  each  u  €  Q\  do 

Qi-QWi «}; 

for  each  (u,v)  €  FS(u)  such  that  d'  +  /«»  <  dj  do 
a  «—  dJ(mod  /mar  +  1); 
dj  4 —  d*  +  luti  > 
fc  «—  dJ(mod  Imax  +  1); 

Pt  — 

Qi  —  <5S\{v}; 

Ql  -  u  {v>; 

end 

7”  -  7"  U  {«}; 
end 

let  z  be  the  next  index  such  that  Q*  ^  0  ; 
for  each  v  £  Q\  do 

for  each  (u,  v)  €  BS(v )  such  that  dj,  +  tuv  <  dj,  do 
a  <—  dJ,(mod  Imax  +  1); 

*—  dj,  +  /«« ; 

6  ♦-  d5,(mod  /mar  +  1); 


K 


»; 


end 

T*  —  T*U{v}; 
end 

endwhile 

comment:  stopping  criterion  met 
0  —  min{d*  +  rf|  :  i  6  T'  U  T*}; 
7-{i€7'JU7T‘  :d?+dj=/?}; 
end 


1.6  The  two-tree  S2  algorithm 

As  in  the  previous  two-tree  algorithms,  the  two-tree  S2  algorithm  uses  mirror  S2  data  structures  to  build 
two  shortest-path  trees.  At  first  one  would  expect  the  stopping  procedure  to  be  the  same  as  in  the  previous 
two-tree  algorithms,  namely,  when  a  node  is  in  both  trees,  find  the  minimum  d\  +  d|  for  all  i  €  T*  U  T*. 
However,  at  the  time  a  node  is  first  placed  in  its  second  tree,  we  are  not  quite  ready  to  search  for  such  a 
minimum  doubly  labelled  node.  In  Nicholson  (1966)  such  a  node  is  proven  to  be  on  a  shortest-path  because 
each  node  in  both  trees  has  had  its  arc-list  fully  scanned.  In  the  S2  implementation  for  each  tree  this  is 
not  the  case.  In  fact,  this  is  the  advantage  of  the  one-tree  S2  algorithm.  In  two  directions  we  must  perform 
additional  scanning  to  meet  the  Nicholson  criterion,  however,  we  will  not  need  to  manage  the  linked-lists. 
All  that  is  needed  is  to  update  distance  labels  and  predecessors.  Actually  we  need  only  scan  a  subset  of  arcs 
from  each  arc-list.  Nicholson  proves  that  any  node  that  is  not  in  either  tree  when  a  node  is  first  in  both  trees 
will  not  be  on  a  shortest-path.  If  arcs  were  scanned  to  these  nodes,  they  would  still  not  be  in  either  tree  or 
on  a  shortest-path,  so  updating  their  distance  labels  would  be  wasted  effort.  These  nodes  also  make  up  the 
majority  of  the  arc-lists.  The  only  arcs  that  need  to  be  scanned  during  the  “mop-up”  phase  are  those  arcs 
that  have  from  nodes  in  one  tree  and  to  nodes  in  the  other  tree.  Since  these  arcs  may  be  scanned  from  either 
node,  it  is  more  efficient  to  consider  these  arcs  from  the  nodes  in  the  smaller  tree.  After  updating  distance 
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labels  and  predecessors  we  are  ready  to  search  for  the  minimum  doubly  labelled  node  to  find  our  minimum 
distance  from  s  to  t  and  a  shortest-path.  The  algorithm  may  be  stated  as  follows: 

Procedure  S22(s,t) 
begin 

initialize: 

pj  -  0 ,df  4-  oo  ,/*(«)  -  1  ,//>(.)  =  |F5(«)|  for  aU  »  €  N ; 

Q‘  -  0  for  z  =  1, ....  I max;  -  {s},  d\  -  O.pJ  *—  s,  T*  «—  0; 

pj  4-  O.dJ  -  oo,  6i(i)  «-  l,6/»(i)  =  |B5(t)l  for  aU  t  €  N; 

Q\  —  0  for  z  =  1, ....  Imax;  Q'a  {t>,  —  O.pJ  —  i,  7*  <—  0; 

while  r*  nr*  =  0  do 

let  z  be  the  next  index  such  that  Q*  0  ; 
for  each  u  6  do 

comment:  determine  next  node  in  FS(u) 

Q\  - 

comment:  determine  new  candidate  arc 
SINSERT(u); 

T*  4-r*  U{u}; 
if  v  ^  O’  then  SINSERT(u); 
end 

let  z  be  the  next  index  such  that  Q\  ^  0  ; 
for  each  v  €  QJ  do 

comment:  determine  next  node  in  BS(u) 

u  «—  wbt(vy, 

comment:  determine  new  candidate  arc 
TINSERT(v); 

T*-V  U{u}; 

if  u  £  Q‘  then  TINSERT(u); 
end 

endwhile 

comment:  mop-up  phase  from  smaller  tree 
if  |T'|<  jT*|  then 
for  each  u  £  T‘  do 

for  each  (u,y)  €  FS(u)  do 
if  y  £  T*  then 

if  d*  +  /uv  <  dj  then 
dj  •  dy  +  /uy; 

pi  *—  «; 

endif 

endif 

end 

end 

else 

for  each  v  £  T*  do 

for  each  £  BS{v)  do 
if  w  £  T‘  then 

if  <4  +/„,„<  dl  then 
dju  *  dj  -I-  Itw  i 

Pw  — 
endif 
endif 
end 
end 
endif 


comment:  stopping  criterion  met 
0  <—  min{d*  +  d\  :  i  €  T‘  U T*}; 
/-{ier*ur  :rf<  + <*!=/?}; 

end 

Procedure  SINSERT(x) 
begin 

fk(x)  —  fk(x)  +  1; 
y  —  W;jt(r); 


while  /Jb(x)  <  //»(x)  and  cf£  +  /ry  >  dj  do 

/*(*)  —  /*(*)  +  i; 
y  *— 

endwhile 

if  fk(x)  <  fh(x)  then 
pi  —  2-; 
dy  <—  d*  4-  Iry ; 

a  <—  dJ(mod  /max  +  1); 


Ql'-Ql  U{x}; 

endif 

end 

Procedure  TINSERT(x) 
begin 

6fc(x)  <—  bk(x)  +  1; 
y  u'kfcOr); 


while  bk(x)  <  bh(x)  and  d\.  +  lyx  >  cPy  do 
bk(x)  *—  bk(x )  +  1; 

y  wbt(x)> 

endwhile 

if  6i(x)  <  bh(x)  then 
P\  *~x\ 


a  «—  d^(mod  /max  +  1); 


endif 

end 


1.7  The  parallel  two-tree  Dijkstra  algorithm 

It  is  readily  observed  that  in  the  two-tree  shortest-path  algorithms,  the  trees  are  independent  of  each 
other.  The  only  requirement  is  to  check  whether  or  not  a  node  is  in  the  opposite  tree.  This  read-only  step 
causes  no  interference  using  multiple  processors.  This  leads  to  the  simplest  asynchronous  parallel  application 
using  two  processors,  one  for  each  tree.  When  one  processor  recognizes  that  a  node  is  in  both  trees,  it  sets  a 
flag  to  tell  the  other  processor  to  find  the  minimum  doubly  labelled  node  in  its  tree  while  it  does  the  same. 
The  minimum  of  the  two  is  the  minimum  path  distance  from  s  to  t.  Again,  a  shortest-path  is  implicit  in 
the  predecessor  labels  beginning  with  the  minimum  doubly  labelled  node.  In  the  algorithms  below,  each 
processor  has  its  own  indentifying  number  called  proad.  The  parallel  processing  construct  called  fork(S) 


indicates  that  two  processors  are  to  be  used  to  execute  the  sections  until  the  join  construct  is  reached.  The 
algorithm  may  be  stated  as  follows: 

Procedure  PD2(s,<) 

begin 

flag  0; 

comment:  begin  use  of  two  processors 
fork(2) 

if  procid  —  1  then 
initialize: 

p‘  *—  Q,d‘  <—  oo  for  all  i  €  N\  Q‘  *—  0; 

synchronize  processors 
while  flag  =  0  do 
for  each  u  £  R‘  do 

for  each  (u,  t;)  S  FS( u)  such  that  d*  +  /« „  <  dj  do 
d*v  * —  d*  +  luv ; 
pi  — 

if  v  ^  Q‘  then  Q‘  <—  Q3  U  {u}; 
end 

T*  <—  TJ  U  {u}; 
if  u  €  T*  then  flag  *—  1; 
end 

comment:  search  Qs  for  minimum  label  nodes  and  place  in  R‘ 
a  ♦-  min{d-  :i  €£?'},/?'  -  {»  :  d?=o},  Q'  -  ; 

endwhile 

/?,  —  min{d*  +  dj  :  i  €  T*}; 

endif 

if  procid  —  2  then 
initialize: 

pj  -  J,  dj  -  oo  for  all  i  €  N\Q'  —  0; 
dj  —  0,p{  *-t,r  -O.tf  -  {<}; 

synchronize  processors 
while  //ai;  =  0  do 
for  each  v  £  R'  do 

for  each  (u,  v)  €  BS(v)  such  that  d[  +  luv  <  dj,  do 

dj,  *—  dj,  +  /uv; 

p‘u  <~v ; 

if  u  6  <p'  then  <p‘  —  (?*  U  {u}; 
end 

r  -r  u{c}; 

if  v  €  7"  then  flag  —  1; 
end 

comment:  search  Q'  for  minimum  label  nodes  and  place  m  R ‘ 
a  *— min{dj  :«€(?'},/?'  >- {i  :  d\=a} 
endwhile 

dt  min{df  +  dj  :  *  €  T*}; 
endif 

comment:  end  use  of  multiple  processors 
join  processors 
0  —  minfd,,/?,}; 

J  -{ieT'ur  :  dj  -f  dj=/?}; 
end 
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1.8  The  parallel  two-tree  Si  algorithm 


The  parallel  two-tree  SI  algorithm  is  similar  to  the  parallel  two-tree  Oijkstra  algorithm.  Each  processor 
is  assigned  one  of  the  nodes,  s  or  t,  and  builds  a  tree  using  the  Si  data  structure.  When  a  node  is  found  to  be 
in  both  trees,  a  flag  is  set  to  tell  both  processors  to  find  the  minimum  doubly  labelied  nooe  in  its  respective 
tree.  The  minimum  of  these  two  gives  the  minimum  distance  path  from  s  to  t  with  a  path  implicit  in  the 
predecessor  labels.  The  algorithm  may  be  stated  as  follows: 

Procedure  PS12(s,<) 

begin 

flag  *—  0; 

comment:  begin  use  of  two  processors 
fork  (2) 

if  procid  =  1  then 
initialize: 

p *  —  0,  d*  <—  oo  for  all  i  6  -V;  Q ’  <—  0  for  z  —  1, Imax; 

Q’o  «-{»}.<*:  ^O’Pl  T‘  -0; 

synchronize  processors 
while  flag  =  0  do 

let  2  be  the  next  index  such  that  Q*  0  ; 
for  each  u  €  Q’  do 

Qt 

for  each  (u,  u)  G  FS( ti)  such  that  dj  -f  luv  <  dj  do 
a  *—  dJ(mod  Imax  +  i); 

dl  *—  d„  +  luv ; 

b  «—  dJ(mod  Imax  +  1); 

pi  *-  «; 

Qt  -  Qt  u  {v}; 

end 

r1 «—  u  {«}; 

if  u  G  T*  then  flag  <—  1; 
end 

endwhile 

/?,  «—  min{d’  +  d|  :  i  G  T*}; 
endif 

if  procid  =  2  then 
initialize: 

pj  <—  0,  d*  «—  oo  for  all  i  G  N :  Q\  —  0  for  z  =  1,  ...,/mai; 

synchronize  processors 
while  flag  —  0  do 

let  r  be  the  next  index  such  that  Q\  ^  0  ; 
for  each  v  G  do 

for  each  (ti,  v)  G  such  that  d[  +  luv  <  d*u  do 

a  «—  d^(mod  /mar  +  1); 

^  "b  i 

6  *-  d^mod  Imax  +  1); 

Pu  —  «; 

Q;-Qlu{u}; 

end 

r  -T'U{r}; 

if  u  G  T*  then  //aj  *-  1; 
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end 

end  while 

Pt  «—  min{d*  +  d\  :  i  €  T*}; 
endif 

comment:  end  use  of  multiple  processors 

join  processors 

P  —  min  {P,,Pt)\ 

J  *—  {«  €  T*  U T*  :d;+dJ=/7}; 

end 


1.9  The  parallel  two-tree  S2  algorithm 

As  in  the  previous  parallel  algorithms,  the  parallel  two-tree  S2  algorithm  assigns  one  processor  to  work 
on  the  tree  rooted  at  s  and  another  processor  to  work  on  the  tree  rooted  at  t.  When  a  node  first  appears  in 
both  trees,  a  flag  is  set  to  initiate  the  mop-up  node  scanning  phase.  As  explained  in  Section  1.6,  we  perform 
the  mop-up  scanning  phase  from  the  smaller  tree  only.  To  perform  this  work  with  two  processors  requires 
each  one  to  share  the  same  data,  namely,  distance  labels.  To  prevent  possible  interference  with  each  other, 
parallel  processing  constructs  called  locks  are  used  so  only  one  processor  at  a  time  may  update  a  distance 
label.  Following  this,  the  processors  are  synchronized  and  the  minimum  doubly  labelled  nodes  are  found  for 
each  tree.  The  minimum  of  these  two  gives  the  minimum  distance  path  from  s  to  i  and  a  shortest-path  is 
implicit  in  the  predecessor  labels.  The  procedures  SINSERT(z)  and  TINSERT(r)  are  as  presented  in  Section 
1.6.  The  algorithm  may  be  stated  as  follows: 

Procedure  PS22(s,t) 
begin 

flag  «-  0; 

comment:  begin  use  of  two  processors 
fork(2) 

if  procid  =  1  then 
initialize: 

p-  «—  0 ,d\  —  oo ,fk(i)  *—  1,  fh(i)  —  |FS(i)|  for  all  i  €  N] 

Q’  —  0  for  z  =  1, ...,  Imax]  <?£  —  {s},  d\  —  0,pj  *-  s,  T*  —  0; 
p]  «— 0,  rfj  ♦—  oo,  bk(i)  *—  1  ,bh(i)  =  |BS(*)|  for  all  i  €  N ; 

Q\  <—  0  for  r  =  1, ...,  lmax\  Q'0  <—  {<},  d]  *—  0,p{  —  l,  T*  <—  0; 
synchronize  processors 
while  flag  =  0  do 

let  z  be  the  next  index  such  that  Q’  ^  0  ; 
for  each  u  €  Q,  do 

comment:  determine  next  node  in  FS(u ) 
v  «-  w/k(u )! 

comment:  determine  new  candidate  arc 
SINSERT(ti); 

T‘  —  T*  U  {u}; 

if  v  i  Q*  then  SINSERT(n); 
if  u  €  T*  then  flag  *—  1; 
end 

endwhile 

endif 


if  procid  =  2  then 
initialize: 

p«  4—  0,  cf|  •—  oo,  bk(i)  -  1,  bh(i)  =  |BS(*)|  for  all  i  6  N\ 
Q,,*~ifoi  z  =  1, Imax ;  Q\>  —  { t },  d\  «—  0,p{  * —  t,  7^  —  0; 
synchronize  processors 
while  flag  —  0  do 

let  z  be  the  next  index  such  that  Q\  £  0  ; 
for  each  v  G  Q*  do 

comment:  determine  next  node  in  BS(u) 
u  <— 

comment:  determine  new  candidate  arc 
TINSERT(v); 

7*  *—  T*  U  {v}; 

if  u  g  Q'  then  TINSERT(u); 
if  v  G  T‘  then  flag  <—  1; 
end 

endwhile 

endif 

comment:  mop-up  phase  from  smaller  tree 
synchronize  processors 
if  JT*|  <  \Tt\  then 

comment:  each  processor  works  on  next  unscanned  node  u  £  T* 
for  each  tigT1  do 

for  each  (u,y)  6  FS(u)  do 
if  y  €  T"  U  V  then 
if  d‘u  +  /uy  <  d*  then 
set  locked 

dy  du  +  lUy  J 

Pi  — 

set  unlocked 
endif 
endif 
end 
end 

else 

for  each  i>  €  T1  do 

comment:  each  processor  works  on  ne^t  unscanned  node  v  €  T* 
for  each  (u>,v)  €E  BS(v)  do 
if  weT'Uf  then 
U  <f[  +  lwv  <  d[,  then 
set  locked 

dyj  1  d^  I  Iwv  » 

Pw  *- 

set  unlocked 

endif 

endif 

end 

end 

endif 

synchronize  processors 

if  procid  —  1,0,  •—  min{d*  +  d}  :  i  €  T*}; 

if  procid  =  2,0,  <—  rrin{dj  -f  dj  :  t  €  T 

join  processors 

0  —  min{/?,,/?t}; 

J  *-  {i  €  T*  UT*  :  dj  +  c/J  — /?} ; 

end 
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2.  COMPUTATIONAL  EXPERIENCE 


All  nine  algorithms  have  been  coded  in  FORTRAN  and  run  on  a  Sequent  Symmetry  S81  using  either 
one  or  two  Intel  80386  processors.  Several  factors  affect  the  performance  of  shortest-path  codes.  First,  the 
number  of  nodes  is  important  for  Dijkstra-type  algorithms,  whose  majority  of  work  is  searching  a  node¬ 
length  array  for  a  minimum  label  node.  Second,  the  average  degree  (|A|/|./V|)  is  important  because  the 
Dijkstra-type  and  Si-type  algorithms  must  scan  entire  forward  (or  backward)  stars  each  iteration,  while 
S2-type  algorithms  typically  scan  only  a  subset.  Finally,  the  cost  range  of  a  network  affects  the  length  and 
sparseness  of  the  Q  array  for  Sl-type  and  S2-type  algorithms,  which  are  searched  each  iteration.  The  cost 
range  and  degree  also  affect  the  number  of  nodes  that  tie  with  a  minimum  distance  label,  thereby  reducing 
the  number  of  searches. 

Four  node  levels  (1000,  2000,  3000,  4000),  ten  average  degree  levels  (5,  10,  15,  20,  25,  50,  75,  100,  125, 
150),  and  three  cost  ranges  (1-100,  1-1000,  1-10000)  were  chosen  as  being  varied  enough  to  demonstrate 
which  factors  were  influencing  performance.  The  total  number  of  random  networks  generated  was  120.  Each 
code  solved  twenty  problems  per  network  using  the  same  randomly  generated  s  and  t  nodes,  yielding  a 
total  of  2400  problems.  Each  data  point  in  Tables  1-3  is  the  sum  of  times  (in  seconds)  to  solve  the  twenty 
problems.  Since  the  S2-type  codes  require  sorted  forward  and  backward  stars,  all  codes  were  given  sorted 
forward  and  backward  stars  to  eliminate  this  as  a  relative  factor  among  them.  It  is  debatable  whether  or 
not  the  sorting  time  should  be  counted  against  the  S2-type  algorithms  since  it  requires  sorted  arc-lists.  Here 
we  simply  assume  that  the  data  is  available  in  pre-sorted  order  and  concede  that  if  it  were  not  available,  the 
S2-type  algorithms  would  not  be  appropriate. 

With  nine  codes  and  120  networks,  many  comparisons  and  observations  can  be  made.  We  highlight  the 
major  points  of  interest.  First,  there  is  some  overlap  with  previous  studies  and  we  wish  to  confirm  previous 
results.  As  in  Dial  et  al.  (1979),  we  find  that  SI  is  better  than  S2  on  only  the  smallest  degree  problems.  As 
the  forward  stars  get  larger,  the  savings  of  scanning  only  a  subset  are  realized.  On  four  high  node  number, 
low  degree  networks,  Mohr  and  Pasche  (1988)  found  that  a  D2-type  algorithm  required  about  38ST  of  the 
time  needed  by  a  D1  algorithm.  Here  we  found  that  the  averages  for  D2  were  65%,  49%,  and  23%  of  that 
for  D1  for  the  cost  ranges  1-100,  1-1000,  and  1-10000  respectively. 

PS22,  the  parallel  two-tree  S2  code,  is  the  overall  winner.  It  had  the  fastest  time  on  108  out  of  the 
120  networks,  while  PS12  was  fastest  on  12  networks.  PS22  was  the  fastest  code  when  the  average  degree 
increased  above  10  on  networks  with  1—100  cost  range  and  when  the  average  degree  increased  above  5  on 
networks  with  1-1000  cost  range.  It  was  always  the  fastest  code  when  the  cost  range  was  1-10000.  PD2 
improved  on  the  networks  with  the  lowest  number  of  nodes,  lowest  cost,  and  highest  degree  —  all  factors 
that  reduce  the  number  of  and  time  for  searches  for  minimum  label  nodes  by  increasing  ties.  Were  the  degree 
increased  to  make  these  networks  much  more  dense,  PD2  might  become  more  competitive. 

Overall,  S22  was  the  fastest  sequential  code,  and  was  even  faster  than  PS12  and  PD2.  It  had  the  fastest 
time  of  the  sequential  codes  on  all  120  networks  In  general,  the  time  required  using  two  Sl-type  trees  is 
about  23%,  26%,  and  36%  of  the  time  using  only  one  SI  tree  on  the  1-100,  1-1000,  and  1-10000  cost  ranges 


respectively.  Similarly,  two  S2-type  trees  require  on  average  23%,  27%,  and  37%  of  the  time  required  by  the 
S2  code  on  1-100,  1-1000,  and  1-10000  cost  range  networks,  respectively. 

Dreyfus  (1969)  comments  that  savings  may  accrue  for  two-tree  algorithms  if  the  stopping  criterion  is 
reached  well  before  N/2  nodes  are  permanently  labelled  in  each  tree.  This  is  definitely  the  case  for  random 
networks.  The  one-tree  algorithms  scan  approximately  50.4%  of  the  nodes  in  the  network  until  t  is  placed 
in  the  tree,  while  the  two-tree  algorithms  scan  only  about  4.7%  of  the  nodes  until  one  is  first  placed  in  both 
trees.  That  is,  two-tree  algorithms  scan  about  9.3%  of  the  nodes  scanned  by  one-tree  algorithms,  resulting 
in  the  above  mentioned  savings  in  time. 

It  should  be  noted  that  we  also  solved  the  above  problems  using  the  efficient  label-correcting  code 
THRESH-X2  (see  Glover,  Klingman,  Phillips,  and  Schneider  (1985)).  The  total  times  in  seconds  for  the  1- 
100,  1-1000,  and  1-10000  cost  range  networks  were  877.9,  963.4,  and  976.0,  respectively.  Since  this  algorithm 
solves  the  one-to-all  shortest-path  problem,  it  was  not  included  in  the  tables.  We  can  see,  however,  that  it 
is  more  efficient  to  use  label-setting  algorithms  for  the  one-to-one  problem  since  they  stop  before  the  entire 
tree  is  built. 

When  using  two  processors,  D2  and  S12  parallelize  nicely  as  the  PD2  and  PS12  codes.  PD2  averages 
a  speed-up  of  1.93  over  D2.  PS12  averages  a  speed-up  of  1.93  over  its  sequential  counterpart,  S12.  In  fact, 
on  some  networks  a  speed-up  over  2.0  is  achieved.  This  is  due  to  ties  for  the  minimum  label  node.  The 
sequential  versions  scan  all  nodes  that  tie  in  one  tree  before  moving  to  the  other  tree.  With  two  processors, 
a  node  that  is  first  scanned  from  a  group  with  ties  could  be  the  one  that  is  placed  next  in  the  opposite  tree 
and  the  remaining  tied  nodes  need  not  be  scanned  as  they  are  in  the  one  processor  version.  Less  work  in 
parallel  results  in  a  speed-up  over  2.0.  If  the  sequential  versions  scanned  a  node  from  each  tree  alternatel- , 
the  speed-up  would  be  less  than  2.0,  but  overall  this  was  slightly  slower  than  scanning  all  nodes  that  tied. 

PS22  averages  only  a  speed-up  of  1.40  over  S22.  This  lower  speed-up  is  due  to  the  additional  cost  of 
using  the  parallel  processing  locks  during  the  relatively  lengthy  mop-up  phase.  That  is,  when  one  processor 
has  locked  a  section  of  code,  the  other  processor  waits  idly  until  the  section  becomes  unlocked  before  it  can 
execute  the  same  section. 


3.  SHORTEST-PATH  HEURISTICS 

There  are  times  when  it  may  be  of  interest  to  quickly  find  “good”  paths  in  a  network.  For  example, 
when  finding  a  starting  solution  to  a  network  flow  problem,  paths  may  be  found  to  send  flow  from  sources 
to  sinks  that  do  not  have  to  be  minimum  paths.  “Good”  paths  at  the  start  may  mean  the  minimum  cost 
flow  will  be  found  more  quickly.  We  find  there  is  a  trade-off  between  the  time  to  find  a  path  and  the  length 
of  the  path  relative  to  a  minimum  path. 

We  have  seen  in  Section  2  that  the  S22  code  is  the  overall  fastest  sequential  code  for  finding  a  shortest- 
path  between  two  nodes.  Recall  that  this  algorithm  requires  a  mop-up  phase  to  scan  all  remaining  unscanned 
arcs  in  the  forward  or  backward  stars  of  the  nodes  in  each  shortest-path  tree  before  a  shortest-path  can  be 
found.  This  mop-up  phase  has  been  found  to  take  from  3%  to  70%  of  the  total  time  for  low  degree  to  high 


degree  networks,  respectively.  However,  before  this  mop-up  phase  we  have  a  node  that  is  in  both  shortest- 
path  trees  and  a  path  from  s  to  t  implicit  in  the  predecessor  labels.  It  may  not  be  an  optimal  path,  but  it 
is  likely  to  be  good  and  is  found  much  quicker  on  higher  degree  networks.  The  heuristic  code  H2  used  the 
path  implied  by  this  first  node  in  each  shortest-path  tree. 

Following  the  mop-up  phase  in  the  optimal  S22  code,  there  is  a  search  over  all  nodes  in  each  tree  for 
the  one  with  the  minimum  sum  of  its  distance  labels.  This  same  search  may  be  done  at  the  end  of  heuristic 
H2,  without  doing  the  mop-up  phase,  to  see  if  there  is  a  better  path  than  the  one  implied  by  the  first  node 
in  both  trees.  Heuristic  H3  is  identical  to  H2,  but  performs  this  additional  step. 

Paths  between  s  and  t  are  known  before  a  node  appears  in  both  shortest-pa‘’  .  :es.  Heuristic  code  HI 
uses  the  path  implied  by  the  node  that  first  has  a  finite  distance  label  from  both  s  and  t. 

Tables  4-6  show  the  times  for  twenty  problems  per  network  on  the  same  networks  used  in  testing  the 
optimal  algorithms.  Also  shown  is  the  percentage  this  time  was  of  the  optimal  S22  algorithm.  As  expected, 
substantial  savings  in  time  are  realized  on  high  degree  networks,  where  the  mop-up  phase  dominates  the  S22 
times.  On  average,  HI  requires  about  65%  of  S22  time,  while  H2  and  H3  require  73%  and  75%  respectively. 
However  on  average,  HI  ranges  from  81%  of  the  S22  time  on  the  lowest  degree  networks  to  46%  of  S22  time 
on  the  highest  degree  networks.  Similarly,  H2  has  an  average  range  of  93%  to  51%  of  S22  time  and  H3  has 
an  average  range  of  95%  to  52%  of  S22  time. 

Tables  7-9  show  how  good  these  paths  are  compared  to  the  actual  shortest-paths.  On  average,  HI  found 
the  shortest-path  55%  of  the  time.  The  average  length  of  its  path  was  8%  greater  than  the  shortest-path 
and  the  worst  path  found  averaged  44%  greater  than  the  shortest-path.  H2  found  72%  of  the  shortest- 
paths  (reaffirming  the  necessity  of  the  mop-up  phase).  Its  average  path  length  was  2.4%  greater  than  the 
shortest-path  and  the  worst  path  it  found  averaged  19%  greater  than  the  optimal  path.  H3  found  93%  of 
the  shortest-paths  with  an  average  path  length  0.5%  greater  than  the  shortest-path.  Its  worst  path  averaged 
5%  greater  than  the  shortest-path. 

It  should  be  noted  that  in  a  few  instances,  HI  found  more  of  the  optimal  paths  for  a  given  network 
than  H2.  (See  Table  10,  nodes  =  1000  and  degree  =  5.)  This  is  similar  to  the  case  in  which  the  first  node 
placed  in  both  shortest-path  trees  is  not  necessarily  on  a  shortest-path.  Here,  the  node  that  first  has  two 
finite  labels  (HI)  is  on  a  shortest-path,  but  is  not  the  node  that  is  first  placed  in  both  trees  (H2). 

4.  CONCLUSION 

The  objective  of  this  paper  has  been  to  present  four  new  shortest-path  algorithms,  two  sequential  and 
two  parallel,  and  to  empirically  compare  them  with  five  algorithms  previously  discussed  in  the  literature. 
The  new  algorithms  combine  the  highly  effective  data  structures  of  the  SI  and  S2  algorithms  with  the  idea 
of  building  vrees  from  a  source  node  and  a  sink  node  in  order  to  find  a  shortest-path.  We  found  that  the 
new  S22  algorithm  was  the  fastest  sequential  algorithm  on  all  networks.  The  new  parallel  algorithm,  PS22, 
was  the  fastest  algorithm  on  all  but  the  lowest  degree  networks,  where  PS12  was  the  fastest.  It  appears  that 
the  parallel  two-tree  Dijkstra  algorithm,  PD2,  might  be  competitive  only  on  very  low  cost,  dense  networks. 


The  secondary  topic  of  this  paper  is  heuristic  S22-type  algorithms  for  obtaining  near-minimum  paths. 
Three  new  heuristic  shortest-path  algorithms  were  discussed  and  were  shown  to  find  very  good  (often  optimal) 
paths  from  a  source  to  a  sink  much  faster  than  the  shortest-path  can  be  found.  These  heuristics  eliminate 
the  time-consuming  mop-up  phase  required  in  the  S22  algorithm  and  are  quite  effective  on  higher  degree 
networks. 


5.  APPENDIX 

Tables  1-3  present  the  computational  results  for  solving  twenty  problems  for  each  of  the  nine  shortest- 
path  algorithms  discussed  above.  Tables  4-6  present  the  computational  results  for  solving  twenty  problems 
for  the  three  S2-type  heuristics  discussed  above.  Tables  7-9  show  how  often  the  heuristics  found  the  optimal 
solution  and  how  far  off  the  solutions  were  when  they  did  not. 


Table  1.  -  -  Time  in  seconds  for  20  problems  (Cost  range:  1-100) 


code 

nodes 

degree 

D1 

SI 

S2 

D2 

S12 

S22 

PD2 

5 

0.79 

6.22 

0.46 

mm 

0.35 

0.41 

10 

1.30 

■fefl 

5.06 

0.61 

Ha 

0.42 

0.43 

15 

1.59 

3.86 

0.69 

wm 

2.18 

0.44 

0.42 

20 

3.65 

2.06 

1.25 

3.64 

0.87 

0.55 

2.07 

0.53 

0.46 

1000 

25 

4.17 

3.06 

1.74 

3.52 

1.01 

0.57 

2.00 

0.61 

0.48 

50 

4.03 

4.57 

2.05 

2.50 

1.50 

0.65 

1.49 

0.82 

0.53 

75 

4.67 

6.21 

2.09 

2.27 

1.91 

0.73 

1.36 

1.05 

0.60 

100 

6.65 

9.00 

4.01 

2.50 

2.63 

0.89 

1.36 

1.30 

0.64 

125 

6.74 

9.87 

4.33 

2.10 

2.57 

0.93 

1.25 

1.36 

0.63 

150 

9.46 

13.72 

5.75 

2.53 

3.32 

1.11 

1.40 

1.70 

0.72 

5 

23.02 

2.28 

2.51 

18.50 

0.92 

0.89 

9.61 

0.59 

0.68 

10 

12.88 

3.23 

2.83 

12.98 

1.13 

0.93 

6.80 

0.69 

0.70 

15 

8.73 

3.43 

2.01 

9.67 

1.34 

0.95 

5.27 

0.78 

0.72 

20 

9.63 

5.50 

3.32 

8.88 

1.55 

0.93 

4.82 

0.83 

0.73 

2000 

25 

9.28 

6.55 

3.35 

8.36 

1.89 

1.09 

4.57 

1.06 

0.78 

50 

10.05 

11.14 

5.31 

5.97 

2.87 

1.21 

3.33 

1.52 

0.92 

75 

9.70 

12.56 

3.65 

4.69 

3.22 

1.23 

2.65 

1.59 

0.89 

100 

13.37 

18.91 

8.79 

4.35 

4.00 

1.41 

2.42 

2.18 

0.98 

125 

mm 

23.21 

9.16 

4.53 

5.18 

1.52 

2.51 

2.45 

1.05 

150 

mm 

27.33 

8.74 

4.81 

5.91 

1.79 

2.62 

2.89 

1.15 

5 

3.23 

3.44 

31.17 

1.30 

15.72 

0.79 

0.92 

10 

5.00 

3.98 

19.70 

1.55 

10.38 

0.90 

0.94 

15 

17.52 

6.86 

4.95 

17.94 

1.96 

1.34 

9.29 

1.09 

1.00 

20 

14.28 

8.14 

5.26 

13.39 

2.22 

1.37 

7.13 

1.21 

0.99 

3000 

25 

14.05 

10.11 

5.96 

12.97 

2.82 

1.47 

6.99 

1.37 

1  08 

50 

11.71 

12.90 

4.56 

8.30 

3.37 

1.53 

4.60 

1.71 

1.11 

75 

14.73 

18.52 

5.32 

7.09 

4.44 

1.73 

4.04 

2.26 

1.23 

100 

20.02 

27.96 

9.09 

7.14 

5.94 

1.96 

3.99 

3.04 

1.28 

125 

22.83 

33.82 

12.85 

6.54 

6.93 

2.14 

3.46 

3.33 

1.33 

150 

23.88 

34.91 

12.52 

5.91 

6.49 

2.14 

3.33 

3.27 

1.32 

5 

48.97 

4.69 

4.83 

42.68 

1.70 

BH 

22.13 

mm 

1.17 

10 

27.81 

6.26 

5.42 

26.38 

2.00 

m 

13.73 

wm 

1.18 

15 

22.66 

8.86 

6.05 

23.14 

2.54 

1.72 

12.77 

1.36 

1.27 

20 

20.28 

11.21 

7.03 

19.05 

2.68 

1.71 

10.11 

1.45 

mm 

4000 

25 

21.14 

15.66 

9.17 

17.76 

3.32 

1.81 

9.39 

1.71 

1.36 

50 

18.55 

20.51 

8.52 

11.65 

4.49 

6.46 

2.18 

iH! 

75 

24.84 

32.71 

14.88 

10.82 

6.69 

5.78 

3.13 

1.55 

100 

27.68 

38.53 

13.72 

9.12 

6.84 

5.03 

3.24 

125 

30.99 

44.23 

13.37 

9.17 

8.66 

2.56 

4.93 

3.90 

1.69 

150 

30.05 

47.07 

■MW 

8.34 

9.63 

2.70 

4.44 

4.39 

1.76 

total 

655.56 

557.49 

120.15 

55.29 

227.58 

Table  2.  -  -  Time  in  seconds  for  20  problems  (Cost  range:  1-1000) 


code 


nodes 

degree 

D1 

SI 

S2 

D2 

S 12 

S22 

PD2 

PS12 

PS22 

5 

30.08 

1.15 

8.70 

0.78 

0.67 

4.62 

0.51 

0.54 

10 

20.47 

1.14 

9.82 

0.87 

0.63 

0.55 

0.53 

15 

20.95 

1.45 

9.37 

0.98 

0.63 

4.97 

0.62 

054 

20 

15.79 

1.38 

8.56 

1.11 

0.63 

4.55 

0.66 

0.55 

1000 

25 

14.87 

3.11 

1.64 

7.60 

1.15 

0.63 

4.07 

0.69 

0.55 

50 

10.59 

4.77 

1.99 

6.91 

1.83 

0.78 

3.76 

1.04 

0.64 

75 

10.74 

7.26 

2.73 

7.04 

2.50 

0.94 

3.71 

1.26 

0.69 

100 

10.54 

9.24 

2.95 

6.12 

2.91 

1.03 

3.33 

1.54 

0.79 

125 

11.78 

11.68 

4.75 

6.22 

3.58 

1.23 

3.36 

1.80 

0.89 

150 

9.90 

10.46 

3.21 

5.07 

3.31 

1.12 

2.76 

1.72 

0.79 

5 

93.53 

2.19 

2.33 

27.70 

1.21 

1.06 

14.14 

0.75 

0.80 

10 

64.37 

3.43 

2.52 

24.63 

1.32 

1.02 

12.65 

0.80 

0.79 

15 

50.13 

4.35 

3.12 

23.35 

1.53 

1.06 

12.18 

0.89 

0.80 

20 

41.88 

5.34 

3.31 

20.52 

1.72 

1.02 

10.69 

0.94 

0.79 

2000 

25 

42.08 

7.38 

4.78 

21.40 

2.07 

1.11 

10.95 

1.08 

0.84 

50 

28.76 

12.52 

5.96 

19.29 

3.35 

1.33 

9.87 

1.70 

0.99 

75 

22.46 

14.42 

6.12 

14.36 

3.85 

1.40 

7.62 

1.90 

1.03 

100 

23.51 

19.21 

4.62 

15.13 

5.39 

1.74 

7.92 

2.56 

1.22 

125 

21.25 

20.12 

5.61 

11.78 

5.02 

1.68 

6.25 

2.55 

1.14 

150 

27.28 

28.97 

8.25 

12.22 

6.68 

2.00 

6.48 

3.13 

1.29 

157.41 

3.44 

3.52 

52.47 

1.67 

1.44 

26.95 

0.99 

■ 

113.51 

5.83 

4.78 

44.69 

1.88 

1.42 

23.37 

1.06 

1 

15 

84.83 

6.93 

4.14 

43.59 

2.14 

1.40 

21.86 

1.18 

20 

65.96 

8.14 

4.81 

42.00 

2.64 

1.52 

21.12 

1  37 

1.13 

25 

63.72 

10.61 

5.97 

39.18 

2.83 

1.50 

19.74 

1.46 

WM 

50 

39.78 

15.73 

7.43 

25.20 

3.84 

1.63 

13.08 

1.88 

Hi 

75 

41.44 

26.88 

11.87 

27.83 

5.98 

2.06 

14.06 

2.85 

HU 

34.30 

26.80 

9.63 

20.02 

5.99 

1.98 

10.48 

2.79 

1.36 

125 

34.88 

32.14 

13.63 

18.56 

6.55 

2.05 

9.72 

3.15 

1.39 

38.17 

39.71 

13.01 

17.74 

7.94 

2.26 

9.30 

3.67 

1.42 

5 

269.68 

4.98 

5.42 

78.45 

2.03 

1.81 

39.54 

1.17 

1.32 

10 

172.27 

8.38 

6.76 

81.07 

2.51 

1.83 

40.44 

1.37 

1.35 

15 

105.57 

7.47 

4.21 

52.22 

2.40 

1.67 

26.15 

1.33 

1.28 

20 

100.07 

11.86 

6.64 

62.87 

3.17 

1.85 

31.88 

1.64 

1.36 

4000 

25 

el. 74 

12.56 

6.36 

56.04 

3.36 

1.89 

28.03 

1.76 

1.41 

50 

54.62 

21.43 

7.01 

40.52 

4.94 

2.05 

20.58 

2.39 

1.48 

75 

51.40 

31.54 

8.64 

36.33 

6.52 

2.36 

18.53 

3.26 

1.65 

100 

51.97 

43.10 

13.61 

34.58 

9.74 

2.93 

17.88 

4.75 

1.81 

125 

39.31 

33.45 

7.67 

26.70 

8.98 

2.84 

13.79 

4.43 

1.84 

150 

oi.74 

50.83 

24.54 

10.33 

2.99 

12.53 

4.82 

1.93 

total 

2223.33 

573.53 

mm 

1090.39 

146.60 

61.19 

558.11 

74.04 

42.72 

Table  3.  -  -  Time  in  seconds  for  20  problems  (Cost  range:  1-10000) 


code 

nodes 

degree 

D1 

SI 

S2 

D2 

S12 

S22 

PD2 

5 

51.21 

3.51 

11.88 

3.86 

2.57 

6.24 

10 

42.60 

2.77 

9.03 

3.01 

1.81 

4.80 

15 

38.74 

2.46 

10.66 

2.97 

1.63 

5.62 

20 

40.55 

3.96 

2.58 

10.70 

2.92 

1.52 

5.62 

1000 

25 

43.22 

4.59 

2.74 

9.99 

2.99 

1.53 

5.30 

50 

30.62 

5.83 

2.53 

9.21 

3.34 

1.42 

4.88 

75 

34.97 

9.09 

4.00 

11.77 

4.31 

1.64 

6.16 

100 

28.98 

10.11 

3.87 

9.32 

4.28 

1.63 

4.93 

125 

24.48 

10.66 

3.02 

9.21 

4.72 

1.70 

4.87 

150 

25.63 

12.82 

3.79 

8.10 

4.79 

1.66 

4.33 

5 

170.47 

5.02 

4.63 

30.39 

4.32 

2.99 

mm 

10 

152.29 

5.32 

4.31 

25.82 

3.57 

2.27 

Wm 

15 

150.98 

6.23 

4.22 

31.87 

3.70 

2.11 

16.38 

20 

129.11 

6.69 

4.15 

33.01 

3.88 

2.03 

16.81 

2000 

25 

109.28 

6.65 

1.19 

26.79 

3.73 

1.90 

13.92 

50 

101.63 

12.15 

5.21 

28.23 

4.87 

2.01 

14.45 

75 

90.39 

17.34 

6.67 

28.80 

5.86 

2.13 

14.80 

100 

66.55 

17.62 

5.15 

27.10 

6.51 

2.27 

14.11 

125 

62.73 

21.42 

5.41 

28.49 

7.54 

2.50 

14.11 

150 

65.88 

28.84 

8.73 

27.35 

8.46 

2.69 

14.11 

5 

380.23 

6.43 

6.07 

53.09 

4.80 

3.44 

27.06 

10 

330.84 

7.48 

5.82 

49.91 

4.15 

2  57 

25.16 

15 

313.76 

9.13 

6.07 

51.80 

4.22 

2.53 

26.47 

20 

309.35 

11.13 

7.24 

50.05 

4.33 

2.37 

3000 

25 

213.96 

5.24 

47.35 

4.42 

2.26 

2-J.08 

50 

165.98 

16.37 

6.01 

52.71 

6.07 

2.43 

26.91 

75 

155.82 

26.33 

13.15 

59.12 

8.16 

2.82 

30.14 

100 

127.05 

29.77 

12.26 

51.31 

8.86 

2.92 

26.22 

125 

81.45 

24.72 

7.66 

35,55 

7.83 

2.54 

17.83 

150 

1^6.22 

41.81 

10.07 

41.62 

9.84 

3.00 

20.97 

5 

741.87 

8.62 

8.65 

108.49 

5.51 

4.02 

10 

386.83 

7.57 

5.51 

61,50 

4.33 

2.92 

KB 

15 

527.01 

12.33 

8.74 

89.12 

4.90 

2.90 

44.47 

20 

416.12 

13.36 

7.33 

82.19 

5.03 

2.75 

41.75 

4000 

25 

381.66 

15.51 

8.32 

67.48 

4.97 

2.64 

33.89 

50 

233.50 

21.68 

6.57 

81.93 

7.14 

2.86 

40.43 

75 

199.77 

29.65 

8.75 

68.66 

8.04 

2.91 

34.36 

100 

187.13 

43.26 

19.67 

69.30 

9.82 

3.19 

34.60 

125 

135.70 

39.82 

11.15 

56.98 

9.93 

3.15 

28.98 

150 

131.56 

47.19 

14.34 

61.04 

12.08 

3.68 

30.88 

total 

6986.12 

617.34 

262.56 

1626.92 

224.06 

98.01 

824.07 

PS12 

PS22 

2.25 

1.97 

1.81 

1.49 

1.77 

1.39 

1.75 

1.32 

1.78 

1.30 

1.95 

1.24 

2.34 

1.37 

2.37 

1.36 

2.57 

1.40 

2.57 

1.39 

2.47 

2.26 

2.11 

1.77 

2.13 

1.68 

2.21 

1.67 

2.15 

1.57 

2.64 

1.60 

3.15 

1.72 

3.43 

1.69 

3.86 

1.89 

4.32 

2.01 

2.54 

2.05 

1.95 

2.47 

1.89 

2.50 

1.86 

3.25 

1.93 

4.29 

2.12 

4.57 

2.16 

4.11 

1.90 

5.04 

2.15 

3.12 

2.94 

2.50 

2.24 

2.78 

2.27 

2.78 

2.15 

2.78 

2.07 

3.81 

2.25 

4.26 

2.24 

5.01 

2.38 

5.16 

2.30 

6.03 

2.60 

I 


Table  4.  -  -  Heuristic  times  in  seconds  for  20  problems  (Cost  range:  1-100) 


degree 

overall 

code 

nodes 

5 

10 

15 

20 

25 

50 

75 

100 

125 

150 

avg 

avg 

1000 

time 

0.33 

0.35 

0.34 

0.34 

0.35 

0.34 

0.34 

23 

0.34 

0.34 

%  of  S22 

71.2 

63.8 

70.1 

62.2 

60.9 

52.3 

47.5 

36.7 

53.7 

2000 

time 

rngji 

0.66 

Pisa 

0.67 

0.67 

HI 

%  of  S22 

cm 

72.9 

70.5 

71.6 

62.6 

56.3 

m 

47.7 

44.5 

37.5 

0.83 

3000 

time 

1*21 

IMS 

B1EE1 

SI'S 

0.98 

HU 

61.7 

%  of  S22 

Bali 

74.3 

m 

64.0 

56.5 

51.1 

45.9 

66.7 

4000 

time 

1.28 

1.31 

iH<ij 

1.29 

1.30 

1.31 

nn 

uu 

1.31 

%  of  S22 

WM 

mi 

76.7 

65.2 

58.4 

55.3 

suS 

66.8 

1000 

time 

0.40 

0.38 

0.41 

mi 

0.39 

0.39 

rngfi 

0.40 

0.40 

%  of  S22 

86.3 

133 

79.0 

74.8 

60.2 

53.4 

46.2 

42.7 

36.4 

62.7 

2000 

time 

0.81 

0.77 

0.79 

0.76 

[•MS 

0.79 

Bip| 

H2 

%  of  S22 

82.9 

83.8 

81.8 

65.2 

52.9 

42.3 

68.2 

0.93 

3000 

time 

mu 

1.12 

1.14 

1.11 

EK3 

1.09 

1.11 

sa 

ifPEl 

1.11 

69.3 

%  of  S22 

IQ 

89.4 

85.6 

80.9 

Em 

63.3 

56.6 

KnEl 

71.9 

4000 

time 

1.44 

1.47 

1.46 

1.48 

fRfcfl 

1.44 

S3 

1.42 

1.46 

1.45 

%  of  S22 

90.7 

85.7 

85.1 

81.5 

Ell 

64.5 

EDUm 

55.5 

54.1 

74.3 

1000 

time 

0.42 

0.44 

0.40 

0.43 

ebi 

0.41 

sm 

0.42 

0.42 

%  of  S22 

90.6 

80.5 

83.2 

79.4 

m 

56.7 

48.8 

ESS 

37.9 

65.9 

2000 

time 

0.83 

0.81 

0.83 

0.79 

0.84 

0.79 

ni-gi 

H3 

%  of  S22 

93.2 

86.5 

87.3 

85.5 

76.6 

67.2 

61.9 

54.7 

43.8 

0.97 

3000 

time 

1.18 

1.19 

1.15 

1.18 

1.17 

1.12 

DBl 

1.12 

mn 

1.16 

71.9 

%  of  S22 

94.9 

95.1 

83.4 

80.4 

76.3 

65.0 

52.5 

74.7 

4000 

time 

1.52 

1.47 

1.51 

Esa 

ngn 

1.47 

1.47 

1.47 

flEEl 

1.49 

$ggjSS! 

%  of  S22 

96.1 

92.7 

87.7 

W 

mu 

73.9 

61.9 

76.1 

Table  5.-  -  Heuristic  times  in  seconds  for  20  problems  (Cost  range:  1-1000) 


degree 

overall 

code 

nodes 

5 

10 

15 

20 

25 

50 

75 

100 

125 

150 

avg 

avg 

1000 

time 

0.44 

0.43 

0.40 

0.40 

0.42 

0.40 

0.41 

0.40 

Esa 

%  of  S22 

73.7 

69.0 

64.0 

62.6 

53.5 

38.7 

33.2 

35.4 

Ea 

2000 

time 

0.86 

[jjggl 

0.76 

BlgJ 

0.75 

0.72 

HI 

%  of  S22 

78.6 

BSD 

68.3 

43.0 

42.7 

37.4 

0.92 

3000 

time 

mi 

mi 

{sin 

1  Kit'll 

nrnii 

1.03 

hni.-K 

mi 

61.4 

%  of  S22 

SZ3 

69.7 

72.7 

64.9 

51.9 

51.9 

51.5 

45.6 

64.6 

4000 

time 

1.42 

HE1I 

1.41 

twm 

rgi 

1.37 

1.37 

ggg| 

1.40 

%  of  S22 

82.8 

77.6 

76.4 

mu 

Ml 

58.1 

47.8 

48.0 

66.2 

1000 

time 

0.59 

0.53 

0.51 

0.47 

0.48 

0.47 

Esa 

0.49 

%  of  S22 

87.7 

84.8 

Egg§ 

urn 

72.8 

60.9 

50.9 

45.5 

64.1 

2000 

time 

twill 

0.87 

0.85 

0.88 

0.88 

nm 

tiBlii 

H2 

%  of  S22 

92.7 

88.9 

82.2 

83.3 

66.1 

59.7 

50.4 

48.3 

42.6 

69.1 

1.04 

3000 

time 

1.34 

1.26 

1.22 

1.23 

ngg 

1.23 

1.15 

1.14 

iltl 

1.21 

%  of  S22 

93.1 

88.4 

87.2 

80.9 

59.5 

58.1 

55.8 

72.7 

4000 

time 

1.68 

1.65 

IB^I| 

1.58 

EEH 

Em 

ra.ri 

ngi 

1.57 

%  of  S22 

92.8 

221 

85.4 

2 

52.5 

gl| 

1000 

time 

0.61 

0.56 

0.52 

ngjn 

0.48 

0.50 

0.49 

0.52 

0.48 

0.52 

%  of  S22 

91.3 

88.7 

83.1 

76.1 

64.2 

48.0 

42.2 

42.6 

67.0 

2000 

time 

0.97 

0.90 

mm 

RlEIil 

mini 

0.91 

niEfi 

SKIS 

H3 

%  of  S22 

95.5 

95.9 

85.0 

85.8 

221 

67.9 

64.3 

52.4 

R!fl 

44.3 

11 

EHU 

3000 

time 

iron 

1.29 

nga 

[H 

1.26 

mi 

1.25 

1.18 

mm 

1.17 

1.24 

72.3 

%  of  S22 

ms 

Hail 

US 

84.1 

i|g|; 

60.8 

59.3 

51.9 

74.4 

CHUM 

time 

1.71 

1.69 

1.53 

1.62 

1.58 

1.55 

1.56 

1.64 

%  of  S22 

94.5 

92.3 

91.8 

87.4 

83  6 

75.7 

66.1 

56.1 

54.4 

53.6 

75.6 

Table  6.-  -  Heuristic  times  in  seconds  for  20  problems  (Cost  range:  1-10000) 


degree 

overall 

code 

nodes 

5 

10 

15 

20 

25 

50 

75 

100 

125 

150 

avg 

avg 

1000 

time 

2.05 

1.54 

mi 

1.23 

1.21 

m 

1.03 

1.00 

1.00 

RjgH 

%  of  S22 

79.8 

85.2 

HI 

80.6 

79.2 

Bo 

63.3 

59.0 

60.4 

ED 

2000 

time 

2.59 

1.86 

m 

1.61 

1.53 

1.42 

1.36 

1.34 

1.35 

1.33 

HI 

%  of  S22 

86.9 

eh 

79.1 

53.9 

49.6 

ED 

1.75 

3000 

time 

2.97 

1.98 

1.81 

1.63 

1.65 

72.0 

%  of  S22 

86.5 

BW1 

83.3 

71.5 

59.6 

57.8 

54.9 

ED 

4000 

time 

3.32 

2.53 

mm 

2.22 

Em 

1.37 

1.98 

1.95 

E2H 

%  of  S22 

82.7 

86.6 

80.8 

58.1 

62.1 

61.8 

54.8 

oo 

1000 

time 

2.46 

1.72 

1.51 

ma 

1.33 

1.15 

till 

1.07 

1.05 

1.39 

%  of  S22 

95.8 

95.2 

92.6 

ESQ 

86.9 

80.6 

63.3 

63.3 

80. G 

time 

2.13 

1.94 

1.83 

1.69 

1.51 

rgv] 

1.48 

1.46 

H2 

%  of  S22 

97.4 

93.7 

92.0 

89.8 

88.6 

Mi 

71.0 

§10 

59.0 

54.4 

78.8 

1.94 

time 

2.54 

2.30 

2.16 

EliTtl 

1.91 

1.92 

TJ6 

1.72 

1.77 

2.15 

79.8 

%  of  S22 

PTSU 

94.8 

91.2 

91.2 

88.6 

78.7 

68.1 

63.6 

67.8 

58.9 

79.9 

4000 

time 

3.87 

2.76 

2.68 

2.49 

2.36 

mi 

2.19 

2.43 

%  of  S22 

96.3 

94.4 

92.4 

90.4 

89.7 

ED 

ho 

67.0 

66.1 

59.6 

79.9 

1000 

time 

2.48 

1.75 

1.53 

1.35 

1.16 

1.19 

1.13 

1.08 

1.42 

%  of  S22 

96.7 

96.6 

94.2 

pH 

88.4 

82.0 

72.6 

69.0 

64.4 

64.7 

82.1 

2000 

2.94 

2.16 

1.97 

1.86 

1.71 

1.58 

1.59 

n 

1.51 

1.50 

1.83 

H3 

98.6 

94.9 

93.5 

91.4 

89.8 

78.7 

74.6 

ED 

60.2 

55.6 

80.3 

1.98 

time 

gm 

MM 

mi 

1.99 

1.96 

1.89 

1.79 

1.84 

Bf||| 

81.4 

%  of  S22 

ED 

m 

wmm 

91.6 

81.8 

69.6 

64.7 

70.2 

61.2 

Sail 

■  ■ 

time 

2.78 

mi 

isgl 

gm 

2.28 

1.56 

2.18 

mi 

2.22 

mm 

%  of  S22 

ED 

95.3 

ED 

ED 

HO 

79.6 

66.1 

68.3 

Mm 

60.5 

W 

Table  7.-  -  Solution  data  for  20  problems  (Cost  range:  1-100) 


degree 

— 

overall 

• 

code 

nodes 

5 

10 

15 

20 

25 

50 

75 

100 

avg 

avg 

%  opt 

85.0 

550 

55.0 

45.0 

70.0 

65.0 

60.0 

65.0 

70.0 

95.0 

66.5 

1000 

%>opt 

2.6 

11.7 

10.9 

11.8 

6.5 

9.4 

7.8 

4.7 

6.0 

1.2 

7.3 

worst  % 

38.6 

57.9 

48.4 

45.9 

50.0 

42  3 

26.7 

30.0 

37.5 

25.0 

%  opt 

EEfil 

50.0 

50.0 

40.0 

65.0 

600 

80.0 

FfiTil 

2000 

%>opt 

7.0 

10.7 

10.7 

8.4 

2.9 

7.6 

4.2 

5.6 

12.8 

7.9 

58  3 

HI 

worst  % 

30.5 

44.6 

71,0 

50.0 

12.5 

28.6 

26.7 

45  5 

38.4 

7.6 

%  opt 

mm 

65.0 

55  0 

45  0 

55.5 

39.5 

3000 

%>opt 

4.6 

9.1 

9.4 

Bill 

12.6 

9.3 

7.2 

7.1 

5.6 

8.0 

worst  % 

29.5 

94.5 

35.1 

53.5 

46.8 

41.7 

25.0 

44.4 

22.2 

%  opt 

gMSl 

70.0 

35.0 

65.0 

RTjJI 

60.0 

65.0 

65.0 

50.0 

■50.0 

55.0 

4000 

%>opt 

3.7 

9.8 

3.4 

4.8 

4.8 

5.2 

11.0 

7.3 

worst  % 

22.9 

46.6 

21.7 

43.9 

20.0 

25.0 

50.0 

57.1 

37.3 

%  opt 

75.0 

85.0 

85.0 

65.0 

80.0 

95.0 

85.0 

ECTSI 

1000 

%>opt 

3.7 

1.0 

2.2 

1.3 

5.0 

3.1 

0.5 

4.8 

3.0 

2.S 

worst  % 

38.6 

8.2 

26.2 

8.1 

23.5 

26-7 

9.1 

33.3 

28.6 

21.7 

%  opt 

65.0 

80.0 

80.0 

60.0 

75.0 

85.0 

90.0 

85.0 

76. 5 

2000 

%>opt 

2.4 

1.7 

2.7 

1.2 

3.1 

2.7 

1.4 

23 

76.1 

H2 

worst  % 

16.8 

22.1 

13.8 

11.3 

16.7 

17.6 

13.3 

12.5 

25.0 

19.2 

2.6 

%  opt 

60.0 

65.0 

|1 

E^'ll 

60.0 

35181 

73.0 

199 

3000 

%>opt 

2.7 

KM 

E  1 

5.0 

0.7 

2.9 

1.4 

2.8 

2.6 

worst  % 

15.7 

37.5 

6.3 

15.4 

9.1 

25.0 

19.5 

1 

%  opt 

65.0 

55.0 

65.0 

191 

70.0 

90.0 

EEfil 

85.0 

75.0 

%>opt 

4.1 

4.5 

2.6 

V 

2.3 

2.4 

2.8 

2.4 

1.7 

2.5 

worst  % 

35.5 

13.7 

20.0 

EEI 

15.2 

16.7 

11.1 

19.1 

%  opt 

100.0 

100.0 

95.0 

100.0 

100.0 

85.0 

100.0 

96.5 

i 

1000 

%>opt 

0.0 

o.o 

0.1 

0.0 

0.0 

1.2 

0.0 

0.6 

0.3 

1 

worst  % 

0.0 

0.0 

1.8 

0.0 

0.0 

7.7 

0.0 

2.7 

%  cpt 

95.0 

100.0 

90.0 

95.0 

195.0 

90  0 

92.0 

2000 

%>opt 

1.7 

0.5 

0.0 

1.8 

1.1 

0.4 

1.0 

1.2 

0.8 

93.6 

H3 

worst  % 

22.1 

7.2 

0.0 

12.1 

133 

5.6 

25.0 

9.1 

9.4 

0.6 

%  opt 

95.0 

95.0 

90.0 

5 1701 

Sgg|I 

90.0 

95.0 

93.0 

7.8 

3000 

%>opt 

0.4 

0.5 

0.3 

ESI 

:3a 

1.3 

0.6 

0.7 

worst  % 

HI 

6.5 

6.2 

ill 

37,5 

14.3 

9.1 

9.6 

%  opt 

90.0 

EEfil 

|iara 

100.0 

93.0 

4000 

%>opt 

0.9 

0.4 

0.9 

| 

0.0 

1.4 

1.2 

0.7 

worst  % 

11.2 

9.6 

9.6 

4.3 

EEI 

0.0 

18.2 

11.1 

9.6 

B-25 


Table  8.-  -  Solution  data  for  20  problems  (Cost  range:  1-1000) 


degree 

overall  | 

code 

nodes 

5 

10 

15 

20 

25 

50 

75 

100 

125 

150 

avg 

avg 

%  opt 

55.0 

60.0 

35.0 

60.0 

55.0 

50.0 

55.0 

45.0 

50.0 

50.5 

1000 

%>opt 

6.6 

8.1 

12.5 

10.9 

4.8 

11.2 

12.8 

6.3 

12.9 

10.4 

9.7 

worst  % 

25.7 

44.7 

64.4 

165.8 

20.3 

40.7 

129.0 

35.2 

73.7 

70.7 

67.0 

%  opt 

55.0 

70.0 

35.0 

70.0 

50.0 

50.0 

55.0 

nsfii 

65.0 

40.0 

53.5 

2000 

%>opt 

4.9 

5.5 

njjpl 

4.2 

5.3 

11.9 

5.6 

8.4 

11.3 

13.1 

52.6 

HI 

worst  % 

19.5 

74.1 

34.3 

31.8 

37.4 

57.1 

23.4 

38.2 

73.7 

82,8 

47.2 

8.2 

%  opt 

55.0 

70.0 

70.0 

35.0 

65.0 

rrotii 

58.0 

46.9 

3000 

%>opt 

FwEl 

2.7 

5.3 

8.6 

5.4 

8.1 

6.1 

worst  % 

48.4 

14.3 

20.7 

35.4 

m 

23.7 

29.1 

21.6 

46.5 

34.3 

%  opt 

25.0 

45.0 

60.0 

60.0 

50.0 

70.0 

48.5 

4000 

%>opt 

10.9 

11.3 

7.4 

10.4 

7.8 

E 

6.5 

Su 

k  : 

5.8 

8.8 

worst  % 

39.1 

36.4 

68.6 

24.3 

53.0 

26.9 

40.4 

38.4 

41.0 

42.9 

%  opt 

65.0 

75.0 

70.0 

85.0 

70.0 

70.0 

1000 

%>opt 

2.8 

5.7 

4.5 

4.9 

2.3 

1.8 

4.1 

0.8 

4.9 

3.3 

worst  % 

38.0 

32.9 

49.4 

49.8 

31.5 

12.3 

12.8 

8.2 

31  3 

27.4 

%  opt 

70.0 

70.0 

Mi 

!£!<■ 

75.  n 

70  n 

%>opt 

2.0 

1.4 

2.3 

2.4 

3.4 

2.3 

1.1 

2.0 

69.9 

H2 

worst  % 

14.6 

17.4 

13.7 

15.8 

12.5 

22.4 

17.2 

EEu 

;;  | 

m 

14  4 

2.4 

%  opt 

65.0 

65.0 

65.0 

EEB1 

70.0 

80.0 

72.5 

EOT] 

3000 

%>opt 

2.9 

1.8 

tfM 

3.2 

E 

1.4 

1.5 

1.9 

j 

worst  % 

27.3 

15.4 

19.9 

41.7 

W 

13.0 

30.3 

11.5 

9.3 

20.1 

! 

%  opt 

65.0 

45.0 

rfilfli 

60.0 

65.0 

80.0 

Slilil 

70.0 

85.0 

fiinlil 

| 

4000 

%>opt 

1.8 

3.7 

5.4 

2.5 

2.3 

1.8 

0.4 

3.7 

2.9 

0.2 

2.5 

worst  % 

14.9 

25.9 

59.2 

18.9 

4.8 

5.2 

26.4 

21.3 

1.3 

20.7 

%  opt 

100.0 

95.0 

90.0 

HSI 

85.0 

95.0 

95.0 

85.0 

1000 

%>opt 

0.0 

0.1 

0.5 

0.9 

0.2 

0.2 

2.8 

0.5 

worst  % 

0.0 

3.0 

8.4 

13.9 

3.9 

4.6 

31.3 

6.6 

%  opt 

80.0  j 

wm 

95.0 

illfl 

85.0 

IjrTjfjM 

100.0 

90.0 

ma 

2000 

%>opt 

0.2 

liB 

2.5 

0.9 

0.0 

0.6 

0.6 

92.1 

H3 

worst  % 

ESI 

5.2 

4.5 

1.3 

20.7 

1.9 

19.1 

0.0 

7.2 

EEU 

%  opt 

90.0 

95.0 

ffilill 

90.0 

95.0 

90.0 

95.0 

95.0 

MjM 

94.0 

7.6 

3000 

%>opt 

1.6 

1.5 

0.2 

1 

0.1 

■ 

0.5 

worst  % 

27.3 

1.2  | 

19.9 

3.1 

20.3 

km 

1.2 

5.2| 

EEI 

8.3 

%  opt 

80.0 

giin| 

90.0 

95.0 

95.0 

95.0 

95.0 

88.5 

4000 

%>opt 

1 

KIM 

0.2 

M 

HI 

0.0 

0.1 

0.6 

worst  % 

m 

HI 

2.0 

21 

0.9 

21.8 

28.9 

1.3 

8.4 

B-26 


Table  9.-  -  Solution  data  for  20  problems  (Cost  range:  1-10000) 


degree 

overall 

code 

nodes 

5 

10 

15 

20 

25 

50 

75 

100 

125 

150 

avg 

avg 

%  opt 

40.0 

45.0 

55.0 

50.0 

50  0 

45.0 

75.0 

45  0 

50  5 

1000 

%>opt 

8.4 

Ell 

9.8 

7.3 

6.9 

14.3 

11.8 

3.2 

7.7 

011! 

9.7 

worst  % 

32.0 

90.7 

40.6 

35.9 

54  9 

48.1 

26.2 

21.9 

59.4 

67.0 

%  opt 

45.0 

55.0 

40.0 

35.0 

70.0 

65.0 

55.0 

45  0 

53.5 

2000 

%>opt 

5.4 

6.2 

[gP 

138 

11.6 

8.1 

6.2 

11.2 

9.3 

8.0 

52  6 

HI 

worst  % 

20.5 

30.8 

29.6 

Saks 

35.1 

52.3 

55.8 

46.0 

40.0 

18.5 

47.2 

8.2 

%  opt 

50.0 

50.0 

65.0 

45.0 

65.0 

jjpisr 

40.0 

65.0 

50.0 

58.0 

46  9 

3000 

%>opt 

9.1 

Els 

7.5 

9.6 

12.8 

5.7 

12.6 

9.0 

6.9 

6.1 

worst  % 

58.8 

22.2 

52.6 

65.3 

42.6 

Eli 

56.5 

99.6 

39.0 

30.5 

%  opt 

40.0 

g0Ei 

60.0 

60.0 

60.0 

48.5 

4000 

%>opt 

8.5 

7.1 

2.9 

8.3 

lliii 

8.5 

6.5 

5.3 

8.8 

worst  % 

45.4 

35.3 

34.5 

35.4 

40.4 

28.7 

42.9 

%  opt 

60.0 

55.0 

45,0 

80.0 

65.0 

75.0 

75.0 

80.0 

70.0 

) 

1000 

%>opt 

3.3 

2.5 

3.9 

1.9 

2.3 

0.7 

2.3 

1.3 

3.3 

1 

worst  % 

12.3 

13.3 

16.2 

20.8 

27.5 

3.5 

11.6 

11.1 

27.4 

i 

%  opt 

80  0 

75.0 

80.0 

750 

75.0 

65.0 

650 

wilii 

70.0 

%>opt 

1.2 

1.2 

3.0 

4.9 

3.3 

1.5 

1.6 

3.2 

5.6 

2.2 

69.9 

H2 

worst  % 

12.5 

10.0 

17.4 

31.8 

10.9 

15.9 

21.6 

65.9 

23.5 

14.4 

2.4 

%  opt 

75.0 

75.0 

65.0 

65  0 

65.0 

72.5 

20.7 

3000 

%>opt 

2.1 

2.2 

3.5 

1.7 

14 

4.2 

2.5 

1.5 

1.9 

worst  % 

11.7 

20.4 

21.1 

4.4 

11.4 

33.1 

17.6 

14  3 

20.1 

%  opt 

60.0 

75.0 

70.0 

65. 0 

80.0 

70.0 

70.0^ 

67.0 

( 

4000 

%>opt 

1.4 

2.5 

2.3 

2.7 

2.7 

0.4 

2.4 

2.4 

2.5 

! 

worst  % 

5.9 

27.1 

16.3 

15.9 

14.7 

5.2 

19.2 

17.5 

20.7 

j 

I 

%  opt 

85.0 

95.0 

85.0 

95.0 

85.0 

85.0 

100.0 

95-0 

94.0 

1000 

%>opt 

0.8 

0.3 

0.4 

0.2 

0.4 

0.8 

0.0 

0.0 

0.5 

worst  % 

6.1 

5.2 

4.5 

2.7 

2.4 

6.3 

0.0 

0.3 

2.2 

6.6 

%  opt 

75.0 

80.0 

90.0 

95.0 

95.0 

95.0 

92.0 

i 

2000 

%>opt 

1.2 

3.1 

0.4 

M 

0.6 

92.1  j 

H3 

worst  % 

EB 

1.6 

6.4 

18.5 

7.9 

liiE 

2.6 

6.9 

0.6  | 

%  opt 

85.0 

100.0 

100.0 

Bll| 

100.0 

95.0 

ijMlM 

90.0 

7.6 

3000 

%>opt 

0.0 

0.0 

1 

0.0 

0.2 

worst  % 

0.0 

0.0 

0.0 

13.3 

1.7 

8.3 

%  opt 

95.0 

90.0 

95.0 

85.0 

jjtojtnl 

90.0 

900 

95.0 

88.5 

%>opt 

0.1 

0.0 

0.2 

1.3 

'  ■ 

0.9 

1.2 

0.4 

0.6 

worst  % 

2.4 

0.2 

3.4 

15.9 

12.8 

17.5 

9.0 

8.4 
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ABSTARCT.  Minimum  cost  network  flow  problems  with  a  piecewise  linear  convex  cost  function  are  t 
model  various  optimization  problems.  They  are  also  used  extensively  to  approximate  nonlinear  cost  fui 
which  may  otherwise  be  difficult  to  handle.  Solving  the  piecewise  linear  problems  using  a  reform 
approach  is  possible  but  may  be  inefficient.  In  this  paper  we  discuss  a  specialized  direct  approach  ; 
implementation  for  solving  such  problems.  The  direct  approach  handles  the  piecewise  linear  structure 
cost  function  by  allowing  each  arc  to  have  varying  costs  on  different  segments.  Computational  result 
been  reported,  from  which  we  conclude  that  using  such  an  approach  has  a  distinct  advantage  over  u 
reformulation  approach. 


1.  Introduction 

Minimum  cost  network  flow  problems  with  a  convex  piecewise  linear  cost  function  (PLN P)  occur  nat 
when  costs  change  abruptly.  This  may  be  due  to  the  imposition  of  penalties  when  certain  limits  or  < 
are  exceeded  by  the  flows  in  the  arcs  (see  Rockafellar[1984j).  For  example,  penalties  can  be  impos 
understocked  or  oversupplied  goods  or  an  increase  in  transportation  cost  due  to  the  load  exceeding  a  sp' 
limit.  Random  demand  in  a  transshipment  network  may  result  in  such  a  cost  function  (Sun  and  Tsai  [] 
Apart  from  direct  applications,  piecewise  linear  functions  are  often  used  to  approximate  functions  th 
more  difficult  to  handle,  and  therefore,  solution  procedures  for  PLNP  are  of  considerable  interest. 

PLN P  can  be  readily  solved  by  a  reformulation  approach.  In  this  approach  each  piece  or  cost  se; 
for  every  arc  in  the  network  is  replaced  by  an  arc  in  an  equivalent  linear  network  flow  problem  (N  P) 
resulting  N P  can  be  solved  by  any  implementation  of  the  bounded  simplex  method  specialized  for  nr 
flow  problems,  but  it  is  a  much  larger  problem  to  solve.  As  the  number  of  pieces  in  the  PLNP  increa 
does  the  dimension  of  the  N P  resulting  from  the  conversion. 
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The  motivation  of  our  study  is  to  be  able  to  efficiently  solve  piecewise  linear  approximations  of  no 
cost  network  flow  problems.  The  approximation  may  have  to  be  revised  several  times,  and  therefore, 
it  using  <»  reformulation  approach  is  undesirable.  Such  an  approach  is  not  only  likely  to  take  longer 
the  increased  dimensions  but  will  also  result  in  additional  burden  due  to  the  data  manipulation  n cede 
each  revision. 

The  growth  of  dimensions  can  be  avoided  if  a  solution  approach  that  can  solve  the  PLNP  dirt 
its  original  form,  is  used.  Fourer  [1985,1988]  and  Premoli  [1986]  presented  direct  approaches  for  : 
piecewise  linear  programs.  Fourer  [1985,1988]  gives  a  comprehensive  discussion  of  such  an  approach, 
approaches  can  be  specialized  for  problems  with  a  network  structure.  Sun  and  Tsai  [1990]  impier 
Fourer’s  approach  for  network  flow  problems.  The  approach  uses  a  bound  on  the  true  reduced  cost 
nonbasic  arcs.  If  none  of  the  nonbasic  arcs  price  favourably  using  the  bound  then  an  optimal  solu 
guaranteed.  The  drawback  is  that  when  a  nonbasic  arc  does  price  favourably,  its  true  reduced  cost 
be  determined  to  ensure  that  its  entry  into  the  basts  will  result  in  an  improvement  of  the  objective  fui 
This  entails  cycle  traces  in  the  network  before  each  pivot  can  be  made. 

In  this  paper  we  outline  a  slightly  different  approach  for  solving  the  PLNP  directly.  Our  approac; 
not  require  the  determination  of  the  true  reduced  cost  of  a  seemingly  favourable  nonbasic  arc,  and  the 
cycle  traces  are  not  needed.  From  our  empirical  results,  we  conclude  that  this  illows  a  significant  red 
in  the  time  taken  to  solve  the  PLNP .  Section  2  gives  an  introduction  to  the  problem  followed  by  a  d« 
discussion  of  the  direct  method  with  our  approach  in  Section  3.  Our  implementation  is  discussed  in  S 
4  and  the  computational  results  are  reported  in  Section  5.  Piecewise  linear  modification  of  30  bend 
problems  given  by  Klingman  et  al.  [1974]  are  used  for  the  computational  testing  and  their  dimensio 
given  in  Table  III.  The  problems  are  solved  using  both  the  reformulation  and  our  direct  approach.  Sign 
reduction  in  computational  times  is  achieved,  as  reported,  when  the  latter  approach  is  used. 


2.  Problem  Definition 


The  problem  of  interest  can  be  defined  as  follows: 


min  f(x)  = 

3  =  1 


s.t.  Ax  =  b 


(P. 


03 


0  <  z  <  u 

where  A  €  3?mxn  is  a  given  node-arc  incidence  matrix,  6  £  9?m  is  the  given  requirement  vector,  u  f. 
the  given  upper  bound  vector,  and 


{Cj ,  if  0  <  Xj  <  uj ; 

; 

c;\  if 

The  arc  ej  has  Sj  segments.  The  breakpoints  C,  uj, . . . ,  uj1  are  such  that  0  <  vj  <  u? ...  <  u’’  =  u 
capacity  of  the  arc  ej).  The  costs  cj, . . cj’  on  the  Sj  segments  of  arc  are  increasing,  i.e, 

0  <  cj  <  cj  . . .  <  cj’  <  oo 

In  the  reformulation  approach  the  j,h  arc  would  be  replaced  by  a  total  of  s}  arcs.  We  solve  the  pr< 
directly  and,  therefore,  retain  its  current  dimension. 

3.  The  Algorithm 

The  algorithm  can  be  seen  as  an  extension  of  the  upper-bounded  simplex  method  for  network  program 
problems.  The  upper-bounded  simplex  method  allows  the  flow  in  a  nonbasic  arc  to  be  on  its  lower  or  i 
bound.  In  our  case  each  arc  is  allowed  to  have  several  breakpoints  and  the  flow  in  a  nonbasic  arc  is  all 
to  be  at  any  one  of  these  breakpoints. 

Definition:  A  basic  solution  is  a  partitioning  of  z,  where  xB  are  the  basic  variables  which  can 

any  value  within  their  lower  and  upper  bounds,  and  z,v  are  the  nonbasic  variables,  which  must  take  o 
value  of  one  of  the  breakpoints. 

Every  basic  solution  corresponds  to  a  rooted  spanning  tree  Tb  in  the  network  (see  Kennington 
flelgason  [1980]).  The  jth  arc  ej  is  considered  to  be  on  the  i,A  segment  if  the  flow  x}  in  it  is  such 
uj-1  <  Xj  <  uj.  If  the  flow  on  an  arc  is  eoual  to  zero,  the  arc  is  considered  to  be  on  the  first  segment  if 
and  on  segment  zero  otherwise.  The  effective  ist  of  any  basic  arv  is  unambiguously  determined,  as  the 
of  the  segment  it  is  on.  Using  these  costs,  the  dual  vector  t  can  be  obtained  from  xB  =  eg,  where  B  i 
basis  matrix. 
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Definition:  For  any  path.  P, 


P  =  {st,et,S2,e2,.  •  •  +1}, 


linking  the  from-node  Sj  to  the  to-node  sn+i  in  7b,  the  orientation  sequence  Q{P)  is  given  by 


if  e*  =  (st.st+i); 
if  ek  =  (si+i,Si) 


(see  Kennington  and  Helgason  [1980]). 


3.1  Pricing 


The  price  of  a  nonbasic  arc  is  defined  as 


Pn,  =  cn,  ~  *r  +  *r 
Pn,  -  cn,  + 

where  ir F  and  vr  are  the  duals  of  the  from-node  and  the  to-node  respectively,  of  the  arc  e,v,  .  c~s 
the  slopes  on  the  right  and  left,  respectively,  of  the  flow  in  the  arc  The  pricing  operation  can  be  se 
an  extension  of  that  in  the  upper-bounded  simplex  method.  Unless  a  nonbasic  variable  xs,  is  at  its 
or  upper  bound,  it  can  enter  the  basis  if  either  increasing  or  decreasing  its  value  will  improve  the  objc 
function,  i.e.,  p^  <  0  and  xs,  is  increased,  or  >  0  and  xs}  is  decreased.  Thus,  the  iitl  and  the  (i  - 
segments  of  each  nonbasic  arc  are  priced  to  determine  if  either  of  them  prices  favourably.  In  the  case  v 
the  arc  is  at  its  lower  (upper)  bound,  the  flow  in  the  arc  cannot  decrease  (increase)  any  further,  and  then 
only  the  first  (last)  segment  is  priced. 

Provosition  3.1  For  a  nonbasic  arc  on  its  i,h  segment,  it  suffices  to  price  only  the  iu‘  and  the  (i  4 
segments  of  each  nonbasic  arc.  If  these  segments  do  not  price  favourably,  the  segments  1, 2, . . . ,  (t  -  1) 
2), . . . ,  sj  for  the  nonbasic  arc  tj  will  not  price  favourably. 

Proof :  The  result  clearly  follows  from  the  fact  that  the  piecewise  linear  costs  on  each  arc  are  convej 


This  pricing  strategy  can  be  viewed  as  a  restricted  basis  entry  approach  as  opposed  to  the  pricing  < 
nonbasic  arcs  in  the  reformulation  approach.  Note  that  the  price  p,v,  gives  the  true  reduced  cost  of  th 
e,Vj  only  when  none  of  the  basic  arcs  are  at  breakpoints,  or  equivalently,  the  solution  is  nondegenerate. 
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The  true  reduced  cost  of  a  nonbasic  arc  e.y;  is  given  by: 


dN,  ~  CN,  ~  XI  icB)i{B 
ieB 

where  dy.  is  the  reduced  cost  when  the  flow  in  the  arc  is  to  be  increased,  and  dtj  is  the  reduced  co 
the  flow  is  to  be  decreased. 

Proposition  3.8  When  some  of  the  basic  arcs  are  at  breakpoints ,  >  p^  and  dy  <  py , 

nonbasic  arc  e s,  ■ 

Proof  :  Let  As,  —  Yj  ■  Then 

^JVy  =  ^N,  XI  (cs)‘ ~  XI  (cs)‘ 

(v,),=-t  (y,).=i 

dN,  =  CN,  +  XI  (cb)>  ~  XI  (Cb)>- 

(r,),  =  -l  (Yj).s.  1 

If  some  of  the  basic  arcs  did  not  have  a  possibly  underestimated  cost  assigned  to  them,  ps,  as  defined 
would  give  the  net  change  in  cost  per  unit  change  of  flow  in  the  nonbasic  arc  being  priced  (see  [3]) 
the  assigned  costs  can  only  underestimate  the  actual  costs,  ps,  and  ds,  will  differ  in  costs  of  oni\ 
basic  arcs  tj  which  will  have  an  increment  in  flow  due  to  the  incoming  nonbasic  arc.  Let  Q  denote  the 
indices  of  all  such  arcs.  This  implies  that  for  i  G  H, 

=  Pn,  +  XI  (cb)‘  ~  (cs)i 

dN,  ~  Pn,  ~  XI  (cb)'  ~  (cb)> 

Since  the  piecewise  linear  costs  are  convex,  (cjg);  —  (c^J,  >  0  for  all  i.  Therefore, 

dN,  >  Pn, 
dN,  <  Pn, 

If  >  0  and  <  0,  the  true  reduced  cost  of  the  nonbasic  arc  es,  cannot  be  favourable  as  well.  Tliei 
if  py^  >  0  and  py,  <  ^  lor  all  nonbasic  arcs,  the  solution  is  optimal. 


3.2  Ratio  Test 


The  amount  by  which  the  flow  is  allowed  to  change  on  any  arc  is  tfte  difference  in  the  value  of  the  n 
breakpoint  and  the  flow.  That  is,  the  flow  in  a  basic  arc  is  allowed  to  move  up  or  down  only  on  the  ct 
segment.  This  differs  from  the  strategy  used  in  Fourer  [1985,1988]  and  implemented  by  Sun  and  T 
[1990],  where  at  a  given  iteration,  flows  can  change  over  several  segments.  Our  approach  eliminates  the 
to  trace  cycles  on  the  tree  since  in  case 


or 


Pv  <  0,  but  dt  >  0 


>  0,  but  djji  <  0 


we  allow  a  pivot,  which  is  degenerate. 


Our  computational  experience  demonstrates  that  the  direct  approach  used  by  Sun  and  Tasi  [199 
ineffective  since  the  cycle  traces  consume  more  time  than  that  saved  by  keeping  the  dimension  of  the  pro 
unchanged.  The  times  taken  by  our  implementation  of  the  two  approaches  have  been  compared  in  Tai 

3.3  The  Algorithm 

STEP  0: 

Find  an  initial  basic  feasible  solution  x°  =  [xfljx'v]  and  let  Tb  be  the  basis  tree. 

Calculate  ir  using  tS  =  eg.  Initially  c}  =  cj( x°)  unless  x°  =  0.  in  which  case  c;  =  cj. 

STEP  1:  Pricing 
Let 

'Fi  =  {e/Vj  :  pxi  <  0  and  flow  xs,  in  arc  ev,  is  to  be  increased  } 

'Fj  =  {es1  '■  pji}  >  0  and  flow  x/v,  in  arc  es,  is  to  be  decreased.} 

If  *Fj  U  Vj  =  <j>,  terminate  with  [x^jx^]  as  the  optimal  solution. 


.herwise  select  e,  6  $|  and  set 


,  f  +1,  if  e;  6  'Fj; 
~\-l.  ifeje’Fj. 
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:  Ratio  Test 


9 


9 


Set 

if  ei  € 

\  X]  —  u”  ,  if  (j  €  #2, 

where  U+  >  x;-  and  u~  <  Xj  are  the  closest  irralrpomts  to  the  right  and  left  of  x; , 
respectively. 

A?  =  min  [x,  —  u*-1} 

0,(P)=i  '  1 

A3  =  min  U-z.) 

-0,(P)=4l  1 

where 

u*"1  <  I.  < 

Define 

A  =  min{Ai,  A2,  A3} 

Let  the  minimum  be  attained  for  arc  e; 

STEP  3:  Update  Flows 

Set  ij  =  Xj  +  A 6  and  Xj  =  -  A6 Oi{P) 

If  A  =  Ai  go  to  step  I 
STEP  4:  "IVee  and  Dual  Update 

Replace  ei  in  the  tree  by  tj.  Update  B,  N  and  the  duals. 

Return  to  step  1. 

4.  Implementation 

NETFLO  [1980],  which  is  an  efficient  implementation  of  the  upper  bounded  simplex  method  on  the  g 
has  been  modified  to  implement  the  direct  simplex  algorithm  and  the  modified  code  is  called  RFSBAS 
data  structures  are,  therefore,  an  extension  of  those  used  in  NETFLO.  RESBAS  has  additional  arc-1 
arrays  used  to  store  the  breakpoints  and  costs  for  the  Sj  segments  of  each  arc  e; .  Another  arc-length  ar 
used  to  indicate  the  current  segment  that  each  arc  is  on.  Additional  node-length  arrays  are  also  introi 
to  record  the  closest  breakpotr.ts  above  and  below  the  current  flow  ij  for  the  basic  arc  e}. 
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The  strategy  outlined  by  Fourer  in  [1988]  and  implemented  by  Sun  and  Tsai  in  [1990],  where  flot 
arc  is  allowed  to  change  over  more  than  one  segment  has  also  been  implemented  for  a  comparative 
The  implementation  is  called  PWNET. 

5.  Computational  Experience 

Standard  network  flow  problems  generated  by  NETGEN"  [1974]  were  modified  to  produce  problem 
piecewise  linear  cost  functions.  The  capacity  of  each  arc  is  divided  into  n  segments  which  are  as 
increasing  costs.  The  number  of  segments  generated  can  be  varied.  All  the  repented  limes  are  the  av 
of  the  time  taken  over  three  runs.  Table  1  lists  the  times  taken  to  solve  the  modified  NETGEN  pie 
linear  problems  with  eight  segments  for  each  arc.  The  problems  were  solved  using  NETFLO,  RESBA 
PWNET.  All  the  three  codes  are  integer  FORTRAN  codes.  To  solve  the  problems  using  NETFLO 
of  the  eight  segments  on  every  arc  is  replaced  by  an  arc  to  yield  an  equivalent  problem  with  a  linea 
function.  The  dimension  of  the  reformulated  problem  is,  therefore,  eight  times  that  of  the  original  pie 
linear  problem.  The  testing  was  done  on  the  SEQUENT  Symmetry  S81  and  the  times  reported  ar 
clock  times  in  seconds.  The  time  taken  for  the  input  of  the  data  and  the  output  of  the  results  is  not  inc. 

Taf  le  I  illustrates  that  on  the  average,  the  direct  approach  used  in  RESBAS  takes  48  percent  of  tin 
taken  by  the  reformulation  approach.  A  large  part  of  the  savings  in  time  is  accounted  for  by  the  fact  th 
the  average,  RESBAS  takes  34  percent  of  the  total  time  in  the  pricing  operation,  while  NETFLO  tal 
percent  of  the  total  time.  The  direct  approach  used  in  PWNET  was  found  to  take  more  time  than  NE' 
in  most  cases.  On  studying  the  time  taken  by  different  parts  of  PWNET  we  concluded  that  the  cycle 
needed  in  this  approach  turn  out  to  be  expensive  and  result  in  its  poorer  performance. 

Table  II  lists  the  comparative  times  taken  by  NETFLO,  RES3AS,  and  PWNET  to  solve  the  NET 
problems  with  80  percent  of  the  arcs  specified  to  be  capacitated.  The  original  NETGEN  problems 
uncapacitated  arcs  in  most  cases.  An  uncapacitated  arc  has  a  large  upper  bound,  and  under  the  absei 
a  capacity  requirement,  the  breakpoints  are  large  as  well.  As  seen  from  Table  II,  having  a  large  perce 
of  capacitated  arcs  does  not  effect  the  performance  of  either  code  considerably.  RESBAS  takes  55  perci 
the  time  taken  by  NETFLO,  on  the  average,  for  these  problems  and  the  performance  of  PWNET  cont 
to  be  poorer  than  that  of  NETFLO  in  most  cases.  PWNET  is  not  used  for  further  testing. 

Table  III  lists  the  comparative  times  taken  by  NETFLO  and  RESBAS  to  solve  NETGEN  problems 


60  percent  of  the  arcs  having  a  high  cost  and  50  percent  being  capacitated.  The  original  NETGEN  p: 
have  only  30  percent  of  the  arcs  with  high  costs  in  less  than  half  of  the  problems.  The  increased  per 
of  high  cost  arcs  has  nearly  no  effect  on  the  performance  of  RESBAS.  From  Table  Ill  we  see  that  R 
takes  52  percent  of  the  time  taken  by  NETFLO  for  these  problems  on  the  average 

A  few  problems  were  solved  using  NETFLO  and  RESBAS  to  see  the  effect  of  varying  the  nur 
segments  used.  The  problems  were  solved  using  4,  8,  16,  24,  and  32  segments.  The  results  are  tabui 
Table  IV,  and  figures  1-5  illustrate  the  comparative  time  taken  by  the  two  codes.  As  the  number  of  sej 
used  is  increased,  the  relative  performance  of  RESBAS  is  seen  to  improve.  This  is  expected  since  w 
increase  in  the  number  of  segments  the  dimension  of  the  problem  to  be  solved  by  NETFLO  becomes 
which  results  in  RESBAS  having  an  increasing  advantage  in  the  time  needed  to  solve  the  problem  ( 
average  for  the  five  problems,  RESBAS  takes  approximately  71  percent  of  the  time  taken  by  NETFLC 
the  number  of  segments  used  is  4.  With  an  increase  in  the  number  of  segments  used  to  32,  the  adv 
that  RESBAS  has  over  NETFLO  increases  significantly,  and  the  average  time  taken  is  seen  to  be  o 
percent  of  that  taken  by  NETFLO. 

6.  Conclusions 

From  our  computational  experience  with  problems  of  varying  size  and  characteristics,  we  conclude  i 
direct  solution  procedure  for  PLNPhas  a  definite  advantage  over  a  reformulation  approach.  On  the  a' 
RESBAS  requires  only  50  percent  of  the  time  required  by  NETFLO  to  solve  the  problems  tested.  T! 
for  some  problems  the  improvement  was  seen  to  be  less  striking  than  the  average,  RESBAS  does  no 
longer  than  NETFLO  to  solve  any  of  the  problems.  Since  most  of  the  savings  in  the  solution  time 
from  the  pricing  operation,  RESBAS  is  likely  to  have  an  increasing  advantage  as  the  density  of  the  pr< 
being  solved  increases.  The  number  of  arcs  needed  in  the  reformulation  approach  to  replace  an  arc  i 
original  problem  also  plays  an  important  role  in  the  extent  to  which  the  direct  approach  out  perforn 
reformulation  appraoch.  As  seen  from  tables  4-8  and  figures  1-5,  the  improvement  in  the  time  tak 
RESBAS  becomes  increasingly  significant  as  the  number  of  segments  used  on  each  arc  is  increased  fi 
to  32. 
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Table  I  TIME*  TAKEN  FOR  ORIGINAL  PROBLEMS 


PRB 

NUM 

NETFLO 

RESBAS 

TIME 

RATIO* 

PWNET 

TIME 

RATIO* 

OPTIMAL 

VALUE 

1 

2.17 

.82 

.38 

2.02 

.93 

209934 

2 

2.28 

.83 

.37 

1.98 

.87 

180820 

3 

3.77 

1.31 

.35 

2.81 

.75 

162949 

4 

3.16 

1.11 

.35 

2.39 

.76 

124369 

5 

4.58 

1.46 

.32 

2.69 

.59 

130774 

6 

5.65 

2.24 

.40 

5.60 

.99 

203737 

7 

9.27 

3.16 

.34 

7.08 

.76 

161329 

8 

9.55 

3.19 

.33 

6.91 

.72 

161288 

9 

11.35 

3.77 

.33 

7.28 

.64 

141902 

10 

11.51 

3.68 

.32 

7.63 

.66 

190070 

11 

2.92 

1.49 

.51 

4.19 

1.43 

8871025 

12 

3.90 

1.49 

.38 

3.23 

.83 

3682451 

13 

2.47 

1.24 

.50 

3.52 

1.43 

7921016 

14 

3.51 

1.39 

.40 

2.99 

.85 

3580569 

15 

4.25 

2.52 

.59 

6.67 

1.57 

10660812 

16 

6.03 

2.43 

.40 

5.14 

.85 

3145173 

17 

4.01 

2.27 

.57 

4.62 

1.15 

i  fli 

18 

4.58 

1.95 

.43 

2.97 

.65 

2814656 

19 

7.02 

4.18 

.60 

9.37 

1.33 

16522487 

20 

11.73 

6.14 

.52 

14.03 

1.20 

19403203 

21 

5.11 

2.78 

.54 

4.52 

.69 

11955635 

22 

8.20 

3.77 

.46 

6.72 

.82 

14674811 

23 

6.10 

313 

.51 

10.12 

1.66 

14471077 

24 

7.11 

3.82 

.54 

10.84 

1.52 

12417227 

25 

9.17 

4.42 

.48 

13  19 

1.44 

9154491 

26 

14.20 

8.10 

.57 

25.17 

1.77 

16534913 

27 

14.59 

7.88 

.54 

28.85 

1.98 

16363553  j 

28 

267.20 

243.92 

.91 

1466.85 

5.49 

104874927 

29 

197.13 

162.15 

.82 

791.95 

4.02 

37115244 

30 

125.51 

71.32 

.57 

290.27 

2.31 

6462989 

Average  Ratio  .48  1.36 

*Wall  clock  seconds  on  the  Sequent  Symmetry  S81  'With  respect  to  NETFLO 
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Table  II  TIME*  TAKEN  WITH  80%  ARCS  CAPACITATED 


PRB 

TIME 

TIME 

OPTIMAL 

NUM 

NETFLO 

RESBAS 

RATIO* 

PWNET 

RATIO* 

VALUE 

1 

3.10 

1.29 

.42 

3.46 

1.12 

234607 

2 

3.84 

1.52 

.40 

3,44 

.90 

218543 

3 

5.06 

1.87 

.37 

3.99 

.79 

203033 

4 

5.43 

1.73 

.32 

3.90 

.72 

164195 

5 

7.76 

2.54 

.33 

4.14 

.53 

108281 

0 

9.46 

3.73 

.39 

9.04 

.95 

296534 

7 

12.38 

4.46 

.36 

8.88 

.72 

204810 

8 

13.50 

5.14 

.38 

8.63 

.64 

139952 

9 

15.85 

5.53 

.35 

10.59 

.67 

131246 

10 

14.73 

4.93 

.33 

8.55 

.58 

136127 

11 

5.35 

3.37 

.63 

8.28 

1.55 

10889306 

12 

10.61 

5.88 

.55 

9.17 

.86 

5186373 

1  13 

4.16 

2.43 

.58 

4.95 

1.19 

9138062 

14 

6.15 

2.84 

.46 

5.44 

.88 

4021478 

15 

7.06 

4.55 

.65 

10.01 

1  42 

12538308 

16 

11.36 

5.88 

.52 

10.93 

.96 

4259552 

17 

5.16 

3.07 

.59 

6.13 

1.19 

8788133 

18 

6.98 

3.25 

.47 

5.21 

.75 

3214784 

19 

1  7.02 

4.18 

.60 

9.37 

1.33 

16522487 

20 

11.73 

6.14 

.52 

14.03 

1.20 

19403203 

21 

5.11 

2.78 

.54 

4.52 

.89 

11955635 

22 

8.20 

3.77 

.46 

6.72 

.82 

14674811 

23 

20.82 

15.35 

.74 

36.38 

1.75 

21816035 

24 

24.37 

17.12 

.70 

40.92 

1.68 

18923270 

25 

31.21 

23.47 

.75 

63.16 

2.02 

17110810 

26 

50.75 

39.14 

.77 

97.90 

1.93 

25418455 

27 

52.81 

44.33 

.84 

111.86 

2.12 

25106519 

28 

356.22 

326.03 

.92 

1646.93 

4.62 

103360085 

29 

269.22 

235.88 

.88 

1024.10 

3.80 

38642574 

30 

158.95 

97.92 

.62 

356.02 

2.24 

7154302 

Average  Ratio 

.55 

1.36 

1 

*Wall  clock  seconds  on  the  Sequent  Symmetry  S81  ’With  respect  to  NETFLO 


Table  III  TIME*  TAKEN  WITH  60%  HIGH  COoT  ARCS 


PRB 

iNUMj 


#  OF 

nodes! 


#  OF 
ARCS 


NETFLO 


RESBAS 


TIME 

RATIO 


OPTIMAL 


1 

200 

1300 

3.04 

1.28 

.42 

274313 

2 

200 

1500 

3.34 

1.33 

.40 

184679 

3 

200 

2000 

4.28 

1.67 

.39 

213952 

4 

200 

2200 

4.53 

1.67 

.37 

162210 

5 

200 

2900 

5.61 

1.90 

.34 

184942 

6 

300 

3150 

9.42 

3.89 

.41 

348180 

7 

300 

4500 

12.42 

4.46 

.36 

245771 

8 

300 

5155 

10.67 

3.82 

.36 

176457 

9 

300 

6075 

16.18 

5.27 

.33 

157974 

10 

300 

6300 

13.83 

4.57 

.33 

199056 

11 

400 

1306 

3.80 

2.18 

.57 

10513520 

12 

400 

2443 

5.44 

2.42 

.44 

3787240 

13 

400 

1306 

3.59 

2.04 

.57 

9568307 

14 

400 

2443 

4.68 

2.02 

.43 

3534618 

15 

400 

1416 

4.26 

2.45 

.57 

9032347 

16 

400 

2836 

7.77 

3.40 

.44 

4211578 

17 

400 

1416 

3.89 

2.05 

.53 

8445932 

18 

400 

2836 

6.33 

2.69 

.42 

3744264 

19 

400 

1382 

7.72 

4.84 

.63 

30869990 

20 

400 

2676 

11.93 

6.22 

.52 

13107867 

21 

400 

1382 

5.30 

2.96 

.56 

19048864 

22 

400 

2676 

8.13 

3.68 

.45 

10090412 

23 

1000 

2900 

15.58 

10.61 

.68 

27878120 

24 

1000 

3400 

13.99 

9.39 

.67 

19301285 

25 

1000 

4400 

23.86 

15.49 

.65 

16501324 

26 

1500 

5107 

34.69 

24.92 

.72 

26560879 

27 

1500 

5730 

37.59 

28.04 

.75 

25442763 

28 

8000 

15000 

317.98 

294.26 

.93 

158984551 

29 

5000 

23000 

218.67 

183.66 

.84 

46185509 

30 

3000 

35000 

135.21 

82.97 

.61 

6713147 

Average  Ratio 


.52 


*Wall  clock  seconds  on  the  Sequent  Symmetry  S8I 
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Table  IV  INCREASING  #  OF  SEGMENTS 


BIliei 

NETFLO 

RESBAS 

TIME 

TIME* 

TIME* 

RATIO 

i 

4 

1.52 

1.01 

.67 

8 

2.18 

.81 

.37 

16 

3.84 

.88 

.23 

24 

5.32 

.87 

.16 

32 

6.52 

.79 

.12 

10 

4 

6.50 

3.75 

.58 

8 

11.51 

3.69 

.32 

16 

19.66 

3.44 

.18 

24 

25.33 

2.76 

.11 

32 

25.67 

2.19 

.09 

15 

4 

2.41 

1.82 

.76 

8 

4.26 

2.51 

.59 

16 

10.28 

4.17 

.41 

24 

13.40 

4.36 

.33 

32 

17.54 

4.87 

.28 

20 

4 

5.30 

3.80 

■7A 
.  *  *• 

8 

6.12 

.52 

16 

10.55 

.35 

24 

60.18 

17.12 

.28 

32 

94.71 

24.32 

.26 

27 

4 

11.96 

9.69 

.81 

8 

14.43 

7.88 

.55 

16 

26.39 

9.93 

38 

24 

41.05 

11.31 

.28 

32 

65.59 

15.73 

.24 

^Wail  clock  seconds  on  the  Sequent  Symmetry  SSI 


Graphs  of  the  Times  Taken  with  Increasing  Number  of  Segments 


Fig  1  SOLUTION  TIMES (#1) 
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